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qq ■ ABSTRACT 

O ■ 

The Planck satellite will observe the full sky at nine frequencies from 30 to 857 GHz. Temperature and polarization frequency 
maps made from these observations are prime deliverables of the Planck mission. The goal of this paper is to examine the effects 
of four realistic instrument systematics in the 30 GHz frequency maps: non-axially-symmetric beams, sample integration, sorption 
cooler noise, and pointing errors. We simulated one year long observations of four 30 GHz detectors. The simulated timestreams 
contained cosmic microwave background (CMB) signal, foreground components (both galactic and extra-galactic), instrument noise 
(correlated and white), and the four instrument systematic effects. We made maps from the timelines and examined the magnitudes 
of the systematics effects in the maps and their angular power spectra. We also compared the maps of different mapmaking codes to 
see how they performed. We used five mapmaking codes (two destripers and three optimal codes). None of our mapmaking codes 
makes an attempt to deconvolve the beam from its output map. Therefore all our maps had similar smoothing due to beams and 
sample integration. This is a complicated smoothing, because every map pixel has its own effective beam. Temperature to polarization 
cross-coupling due to beam mismatch causes a detectable bias in the TE spectrum of the CMB map. The effects of cooler noise 
and pointing errors did not appear to be major concerns for the 30 GHz channel. The only essential difference found so far between 
mapmaking codes that affects accuracy (in terms of residual root-mean-square) is baseline length. All optimal codes give essentially 
indistinguishable results. A destriper gives the same result as the optimal codes when the baseline is set short enough (Madam). For 
longer baselines destripers (Springtide and Madam) require less computing resources but deliver a noisier map. 

Key words. Cosmology: cosmic microwave background - Methods: data analysis - Cosmology:observations 



19 

20 



* The main part of the work reported in this paper was done in May 
2006 when the CTP Working Group of the Planck Consortia met in 
Trieste. The author list reflects the CTP membership at the time. Since 
the Trieste meeting, the CTP group has received new members who 
were not involved in this study and whose names do not therefore ap- 
pear in the author list. 



1. Introduction 

Starting in 2003, Planck Working Group 3 (the "CTP" group) 
undertook a comparison of mapmaking codes in increasingly 
realistic situations. The approach to realism proceeded in four 
steps, named after the locations of working meetings of the 
group (Cambridge, Helsinki, Paris, and Trieste). Results from 
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the Cambridge, Helsinki, and Paris steps have been presented 
in previous papers (Poutanen et al. 2006, Ashdown et al. 2007a, 
2007b ). Here we present results from the Trieste simulations, de- 
signed to determine how mapmaking codes handled four aspects 
of real Planck data not included in previous simulations. The 
first was non-axially-symmetric beams. In previous simulations, 
we assumed that the beams on the sky were axially symmetric 
Gaussians. The second was the effect of detector sample inte- 
gration, which introduces an effective smearing of the sky sig- 
nal along the scanning direction. The third was "cooler noise", 
representing the effect of temperature fluctuations induced in the 
focal plane by the 20 K sorption cooler. The fourth was the point- 
ing errors. In previous simulations, we assumed that the detector 
pointings were known without error in the mapmaking. In this 
paper we present the results of this latest round of simulations, 
which are realistic enough to allow us to draw some preliminary 
conclusions about mapmaking. We also outline additional work 
that must be done before final conclusions can be drawn. 

The organization of the paper is as follows. In Section 2 
we describe the simulations that produced the time-ordered data 
(TOD) streams that were inputs to our mapmaking. In Section 3 
we give the inputs that were used in these simulations. In Section 
4 we describe the mapmaking codes we used in this study. 
Section 4 details the changes that we made in those codes since 
our earlier Paris simulation round. Section 5 gives the results 
of our Trieste simulation round and the computational resource 
requirements of our mapmaking codes are listed in Section 6. 
Finally we give our conclusions and proposal for future map- 
making tests in Section 7. In Appendix A we describe an analytic 
model that we used in explaining the effects of beam mismatch 
in the CMB maps. Appendix A also shows how we can use this 
model to correct these effects from the observed spectrum. 

2. Simulations 

We used the Level-S simulations pipeline (Reinecke at al. 2006 ) 
to generate 1-year intervals of simulated detector observations 
(time-ordered data streams, or TODs). As in the Paris round 
(Ashdown et al. 2007b) all simulations were done at 30 GHz, 
the lowest Planck frequency. This was chosen, because the TOD 
and maps are the smallest in data size for this frequency, mini- 
mizing the computer resources required for the simulations, and 
because the beams are furthest from circular, emphasizing one 
of the effects we are trying to study. We simulated the relevant 
sky emissions (CMB, dipole, diffuse galactic foreground emis- 
sions, and the strongest extragalactic point sources) in both tem- 
perature and polarization, plus a number of instrumental effects: 
uncorrelated (white) noise, correlated (1//) noise, noise from 
sorption cooler temperature fluctuations, both circular and ellip- 
tical detector beams, sample integration, semi-realistic nutation 
of the satellite spin axis, and fluctuations of the satellite spin rate. 

TODs 366 days long were generated for the four 30 GHz 
LFI detectors (Low Frequency Instrument), with 1.028 x 10 9 
samples per detector corresponding to a sampling frequency of 
/ s = 32.5 Hz. 

For every sky component we made four different simulated 
TODs. A TOD included the effects of either axially symmetric 
or asymmetric Gaussian beams and the sample integration was 
either on or off. We call these four TODs as 

- Symmetric beams & no sampling 

- Symmetric beams & sampling 

- Asymmetric beams & no sampling 

- Asymmetric beams & sampling. 



For the instrument noise (uncorrelated + correlated) we used the 
noise TODs of the Paris round (Ashdown et al. 2007b). Finally, 
we had a TOD of the cooler noise. Maps were later made from 
different combinations of these TODs. 



3. Inputs 

3.1. Scanning strategy 

The correspondence between the sample sequence of the TOD 
and locations on the sky is determined by the scan strategy. The 
Planck satellite will orbit the second Lagrangian point (L2) of 
the Earth-Sun system (Dupac & Tauber 2005), where it will stay 
near the ecliptic plane and the Sun-Earth line. 

Planck will spin at ~ 1 rpm on an axis pointed near the 
Sun-Earth line. The angle between the spin axis and the optical 
axis of the telescope (telescope line-of-sight, see Fig. [2j is 85°; 
the detectors will scan nearly great circles on the sky. The spin 
axis follows a circular path around the anti-Sun direction with 
a period of six months; the angle between the spin axis and the 
anti-Sun direction is 7?5. The spin axis thus follows a cycloidal 
path across the sky, (like the one we used in our Paris round, 
Ashdown et al. (2007b)). In this simulation the spin axis is re- 
pointed hourly. During each repointing the projection of the spin 
axis onto the ecliptic moves by a fixed 2.'5. Our simulation had 
8784 repointings in total. We assumed non-ideal satellite mo- 
tion, with spin axis nutation and variations in the satellite spin 
rate. 

The scan strategy planned for flight differs from the one used 
here only in that instead of repointing once per hour with a 
fixed offset in ecliptic longitude, we will repoint in 2' intervals 
along the cycloid. To maintain the 2.'5hr~ 1 average rate of mo- 
tion along the ecliptic, the time spent at a given spin axis pointing 
will vary somewhat. 

Spin rate variations were chosen randomly at every repoint- 
ing from a truncated Gaussian probability distribution with pa- 
rameters (0?1 s~ l RMS, 0?3 s~ l max). The abbreviation "RMS" 
refers to the root-mean-square. 

The satellite spin axis nutated continuously according to the 
satellite dynamics. The nutation amplitude was chosen randomly 
at every repointing to mimic the disturbance that the repointing 
maneuver causes in the spin axis motion. In this simulation the 
Gaussian distribution of nutation amplitudes ranged from 0.'006 
to 16.4, with mean and standard deviation values of 1.'2 and 0.'8. 
The nutation amplitudes of the repointings were all < 3. '2 except 
for two large excursions (out of 8784) of 43 and 16.4. This level 
of nutation is many times larger than is (now) expected in flight. 

The discussion of the performance of the Planck pointing 
system is outside the scope of this paper. Our simulation of 
pointing is based on the results of a detailed simulation of the 
spacecraft pointing dynamics. In the simulations of this paper 
the satellite attitude ("satellite pointing") was sampled at 1 Hz. 
In flight the Planck attitude will be sampled at a higher rate. 

We used the HEALPi)(Q pixelisation scheme (Gorski et 
al. |2005a) with N s ide = 512. A map of the full sky contains 
12/V s 2 ide pixels. The Stokes parameters Q and U at a point on 
the sky are defined in a reference coordinate system (eg, e^, n), 
where the unit vector eg is along the increasing 6 direction, e v is 
along the increasing ip direction, and n points to the sky (Gorski 
et al. I2005bl ). The angles 8 and <p are the polar and azimuth an- 
gles of the spherical polar coordinate system used for the celes- 
tial sphere. 

1 http://healpix.jpl.nasa.gov 
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Fig. 1. Number of hits per pixel («hit) for the scan strategy applied in this study. The hit map is shown in the ecliptic (left) and 
galactic (right) coordinates. The latter map shows the areas of the ecliptic poles more clearly. Both maps include the hits of all four 
LFI 30 GHz detectors. The scale is log 10 («hit)- 



The number of hits per pixel from all detectors is shown in 
Fig- IH At this resolution every pixel was hit. 

3.2. Telescope beams 

The horns of the LFI detectors sit in the Planck telescope fo- 
calplane (Fig.|2]l. The center of the field of view, which is empty 
in the figure, is populated with the beams of the HFI bolometers 
(High Frequency Instrument). There are two 30 GHz horns in the 
focalplane. The corresponding beams are labelled with "27" and 
"28" in Fig. |2] Behind each horn we have two detectors tuned to 
orthogonal linear polarizations, called LFI-27a, LFI-27b, LFI- 
28a, and LFI-28b. For simplicity, we refer hereafter the field of 
view as focalplane. 

The time-ordered data were simulated using two sets of 
beams. The first set were circular Gaussian beams of the same 
beamwidth for all of the detectors. The second set were the best- 
fit elliptical beams for the LFI 30 GHz detectors. The beam pa- 
rameters that we used in our simulations, are given in Table 1. 
These were obtained by fitting a bivariate Gaussian to the co- 
polar component of each beam over the whole angular area in 
which each beam was calculated. For the 30 GHz beams this 
was -0.026 < u,v < 0.02ffl 

Realistic main beams have been simulated in the co- and x- 
polar basis according to the Ludwig's third definition (Ludwig 
119731 1 in Mv-spherical grids with 301 x 301 points (Am = Av - 
10~ 4 ). Each main beam (Fig. [3) has been computed in its own 
coordinate system in which the power peak falls in the center of 
the Mv-grid and the major axis of the polarization ellipse is along 
the M-axis. In this condition, a well defined minimum appears in 
the x-polar component in correspondence to the maximum of the 
co-polar component. 

Main beam simulations have been performed using the phys- 
ical optics considering the design telescope geometry and nom- 
inal horn location and orientation on the focalplane as de- 
scribed in Sandri et al. C2004). The computation was carried out 
with GRASP8, a software developed by TICRA0 (Copenhagen, 
Denmark) for analysing general reflector antennas. The field of 



2 u and v are equal to sin(#) cos(<,p) and sin(#) sin(i/>), where 9 and if 
are the polar and azimuth angles of the spherical polar coordinates of 
the beam coordinate system (see Fig. [2}. The wv-coordinate system is 
applied in the antenna beam pattern representations since it permits to 
map the beam from the spherical surface to a plane. 

3 http://www.ticra.com 
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Fig. 2. Footprint of the Planck LFI focalplane on the sky as seen 
by an observer looking towards the satellite along its optical axis. 
The origin of a right-handed Mv-coordinate system is at the center 
of the focalplane (telescope line-of-sight). The z-axis is along the 
line-of-sight and points towards the observer. Labels "18-23" 
refer to 70 GHz horns, "24-26" refer to 44 GHz horns, and "27- 
28" refer to 30 GHz horns. Each beam has its own coordinate 
system as shown in the figure. The coordinate axes show the 
polarization orientations for the co-polar beams of the "a" and 
"b" radiometers. The focalplane scans the sky as the satellite 
spins. The scanning direction is indicated by an arrow. The +u 
axis points to the spin axis of the satellite. The centers of the 
30 GHz beams sweep -1° from the ecliptic poles when the spin 
axis is in the ecliptic plane. In Appendix [A] we need the angle 
from the scan direction to the M-axis of horn "27". This angle is 
67?5. 



the source (feed horn) has been propagated on the subreflector to 
compute the current distribution on the surface. These currents 
have been used for evaluating the radiated field from the sub re- 
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Fig. 3. Contour plot in the Mv-plane (-0.026 < u,v < 0.026) of the main beam co-polar component computed for the 30 GHz feed 
horns, assuming an ideal telescope. The color scale goes from -90 to dB. The fit bivariate Gaussian contours are superimposed 
with dotted lines. From left to right the beams are for LFI-27a/b and LFI-28a/b. They are perfectly symmetric beams with respect 
to the M-axis because of the symmetry of the Planck LFI optics. 

TABLE 1 

Beams 



Detector FWHM" Ellipticity* if f m ff m ifc. 



Symmetric 

27a 32.'1865 1.0 0?2 ... 153?6074 4?3466 -22?5 

27b 32:1865 1.0 89?9 ... 153?6074 4?3466 -22?5 

28a 32:1865 1.0 -0?2 ... -153?6074 4?3466 22?5 

28b 32:1865 1.0 -89?9 ... -153?6074 4?3466 22?5 



Asymmetric 

27a 32:2352 1.3562 0?2 10L68 153?6074 4?3466 -22?5 

27b 32.'1377 1.3929 89?9 100?89 153?6074 4?3466 -22?5 

28a 32.'2352 1.3562 -0?2 78?32 -153?6074 4?3466 22?5 

28b 32:1377 1.3929 -89?9 79? 11 -153?6074 4?3466 22?5 



a Geometric mean of full width at half maximum (FWHM) of the major and minor axes of 
the beam ellipse. Symmetric beam FWHM was chosen to be the arithmetic mean of the 
two FWHMs of the asymmetric beams. In practice the beamwidths will not be known to 
this level of precision, but we give additional significant figures here to show the level of 
variation of the widths, and to reflect what was actually used in the simulations. 

b Ratio of the FWHMs of major and minor axes. 

c Angle between M-axis and polarization sensitive direction (see Fig.©. 
Angle between t(-axis and beam major axis (see Fig. [2}. This angle is irrelevant for axially 
symmetric beams. 

c Angles giving the position of the detectors in the focalplane. They give the rotation of the 
detector ;/vz-coordinate system from its initial pointing and orientation (aligned with the 
telescope line-of-sight uvz-axes) to its actual pointing and orientation in the focalplane (see 
Fig.©- 



Hector. The calculation of the currents close to the edge of the 
scatterer has been modeled by the physical theory of diffraction. 
The radiated field from the sub reflector has been propagated on 
the main reflector and the current distribution on its surface is 
used to compute the final radiated field in the far field. 

In this paper we considered the effects of the co-polar main 
beams only and did not include the effects of the x-polar beams 
in our simulations. In making maps from data with both circular 
and elliptical beams, we can quantify the effect of the elliptical 
beams in the maps and in the power spectra derived from the 
maps. 

3.3. Signal sampling and integration 

The readout electronics of the LFI 30 GHz channel sample the 
signal measured by the detectors at 32.5 Hz. The value recorded 
in each sample is the average of the measured signal over the 
period since the last sample. This non-zero integration time has 



the effect of widening of the beam along the scan direction. If the 
spin speed remains constant throughout the mission, this effect 
cannot be separated from the shape of the beam. 

To quantify this effect, the TOD have been simulated using 
two options for the sampling. The first is to use instantaneous 
sampling, where the signal is not integrated over the past sample 
period, rather the sample value is given by the sky signal at the 
instant the sample is recorded. This option gives an idealized 
result to compare to the realistic sampling behaviour. 

In the second and realistic option, there is an additional ef- 
fect which must be taken into account. In the Level-S simulation 
pipeline, the pointing of the detector is sampled at the same rate 
as the signal from the detector and given at the instants the sam- 
ples are taken. However, the effect of the integration time is to 
smear the sample over the past sample period; in effect, the re- 
ported pointing lags the signal by half a sample period. In order 
to minimise the residuals in the mapmaking, some of the codes 
used to produce the results in this paper can perform an interpo- 
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lation to shift the pointing back by half of a sample period to the 
middle of the sample. 

3.4. Noise 

3.4.1. Detector noise 

We used the instrument noise from our Paris round of simula- 
tions (Ashdown et al. 2007b). Its uncorrelated (white) noise was 
simulated at the level specified in the detector database. Its nom- 
inal standard deviation per sample time was cr = 1350 /jK (ther- 
modynamic (CMB) scale). Correlated 1 // noise was simulated 
from a power spectrum with a knee frequency of 50mHz and 
slope -1.7. For the details of the noise generation see Ashdown 
et al. I2007bl Subsequent tests of the 30 GHz flight detectors 
show a lower knee frequency than 50mHz, so these simula- 
tions can be taken as providing a conservative upper limit on 
1 // noise. No correlation was assumed between the noise TODs 
of different detectors. For the optimal and Madam mapmaking 
codes, perfect knowledge of the noise parameter values was as- 
sumed in the mapmaking phase. 

3.4.2. Sorption cooler temperature fluctuations 

The Planck sorption cooler has two interfaces with the in- 
struments, LVHX1 with HFI and LVHX2 with LFI (LVHX = 
Liquid- Vapor Heat eXchanger). The nominal temperature of 
LVHX1 is 18K, providing precooling for the HFI 4-K cooler. 
The HFI 4-K cooler in turn cools the HFI housing and the LFI 
reference loads. The temperature of the HFI housing is stabilized 
by a Proportional-Integral-Differential (PID) control. LVHX2 
determines the ambient temperature of the LFI front end. Its 
nominal temperature is 20 K. 

Temperature fluctuations from the coolers affect the LFI data 
in three ways. First, fluctuations in LVHX2 propagate through 
the LFI structure to the LFI horns, resulting in fluctuations in 
additive thermal noise from the throats of the horns (where the 
emissivity is highest). Second, fluctuations in the LFI structure 
driven by LVHX2 propagate to HFI both by radiation and by 
conduction through struts supporting the HFI, and thence the 
LFI reference loads. Any temperature fluctuations of the refer- 
ence loads will appear as spurious signals in the LFI detectors. A 
significant part of these fluctuations are suppressed by the HFI 4- 
K PID control, but the 30 GHz loads are between the struts and 
the control stage, so fluctuations are incompletely suppressed. 
Third, temperature fluctuations of the LFI reference loads are 
also driven by LVHX1, propagated indirectly through HFI. This 
last effect is quite small. 

A coupled LFI/HFI thermal model was not available when 
the work reported in this paper was performed, so we were un- 
able to include all of these effects realistically. Instead, in this 
paper we consider only the direct effect of LVHX2 instabilities 
on the feeds. 

The propagation of LVHX2 temperature fluctuations to LFI 
output signals involves two transfer functions (TF): 

- TF1 describes how the temperature fluctuations of the cold 
end propagate to the temperature fluctuations of the LFI front 
end. 

- TF2 describes how the fluctuations of the ambient temper- 
ature of the LFI front end translate into a variation of the 
output signal of a detector. 

The Planck LFI instrument team developed TF1 from the 
LFI thermal model. TF2 is described in Seiffert et al. ([2002 ). 



The thermal mass of the LFI front end suppresses fast tempera- 
ture variations, so TF1 rolls off steeply at high frequencies. TF2 
is a constant multiplier that is different for detectors of different 
frequency channels. The impact of the sorption cooler tempera- 
ture fluctuations on the output signals of the LFI detectors has 
been discussed by Mennella et al. (2002). 

For this study, simulated cooler TODs for the four LFI 
30 GHz detectors were generated as follows. The LFI instrument 
team applied TF1 and TF2 to a ~~ 100 hour sequence of LVHX2 
cold end temperatures taken during cooler operation, producing 
a -100 hour segment of data that approximates the fluctuations 
as they would appear at the output of an LFI 30 GHz detector. 
Three hours of these data from the beginning and end of the 
chunk are shown in Fig. [4] The cooler signal has a distinct peri- 
odicity, whose cycle time is -760 s. Every sixth peak is stronger 
than the other peaks. Fig. [4] shows that the cycle time and the 
amplitude of the sorption cooler fluctuations remain stable over 
the -100 hour period. 

We used linear interpolation to increase the sampling rate of 
the cooler signal (Fig. |4| from its original 1 Hz to the detector 
sampling rate 32.5 Hz. After that we glued a number of these 
-100 hour segments one after another to obtain a one year long 
cooler TOD. We used a -10 hour overlap in the boundaries of 
the successive segments. The segments were manually adjusted 
in the time axis to give a good alignment of the fluctuation peaks 
and valleys in the overlap region. Finally we multiplied the end 
of a previous segment with linear weights descending from 1 to 
and the beginning of the next segment with linearly ascending 
weights (ascending from to 1) before summing the segments 
in the overlap region. 

The resolution of the LFI thermal transfer function model 
could not distinguish between different detectors at 30 GHz at 
the time of this work. Therefore, in this paper all four LFI 
30 GHz detectors were represented by the same one-year-long 
cooler TOD. 



3.5. Dipole 

The temperature Doppler shift arises from the (constant) motion 
of the solar system relative to the last scattering surface and from 
the satellite motion relative to the Sun. The latter signal is usu- 
ally used for the calibration of the CMB observations. In this pa- 
per we assumed a perfect calibration and did not study the effects 
of calibration errors in our maps. We therefore chose to include 
the temperature Doppler shift of the solar system motion in our 
simulations, but we did not include the part that arises from the 
satellite motion relative to the Sun. 



3.6. CMB 



As in the Paris round (Ashdown et al. 2007b), the CMB template 
used here is WMAP (Wilkinson Microwave Anisotropy Probe) 
constrained as described in the following, and included in the 
version 1 . 1 of the Planck reference sk)Q. It is modelled in terms 
of the spherical harmonic coefficients, a,£' B , where T refers to 



4 The CMB and extragalactic components of the 
Planck reference sky vl.l used here are available at 
http://people.sissa.it/-planck/reference_sky The diffuse Galactic 
components are available at http://www.cesr.fr/-bemard/PSM/ 
The most recent version, named Planck Sky Model, in- 
cluding the CMB template used here, is available at 
http://www.apc.univ-paris7.fr/APC_CS/Recherche/Adamis/PSM/psky-en.html 
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Fig. 4. Variations of the output signal (TOD) of an LFI 30 GHz detector caused by the temperature fluctuations of the sorption cooler 
cold end. We show here the beginning (left) and the end (right) of a -100 hour segment of data. The signal is a result of applying 
LFI thermal transfer functions to the cold end temperature data measured from the sorption cooler flight hardware. The thermal 
transfer functions were derived from the LFI thermal model developed by the Planck LFI instrument team. The vertical axis is 
antenna temperature in microkelvins at 30 GHz. 



temperature, and E and B refer to the polarization modes. The 
a tm were determined for multipoles up to i — 3000. 

For I < 70, the aK were obtained by running the anafast 
code of the HEALPix package on the first-year WMAP CMB 
template obtained by Gibbs sampling the data (Eriksen et 
al . 120041 1. The a E [m were then given by 




(1) 



where Cf Y (X, Y = T, E) is the best fit angular power spectrum 
to the WMAP data, and X( m and yc, n are Gaussian distributed 
random variables with zero mean and unit variance. For m = 0, 
the imaginary part of yc m and the V2 were not applied. 



For 



T,E 



> 70, we used the synfast code to generate the a. 
as a random realization of the Q coefficients of the theoretical 
WMAP best fit cosmological model. 

3. 7. Foreground emission 

With the exception of the Sunyaev-Zel'dovich (SZ) signal from 
clusters of galaxies, foregrounds have been modelled according 
to vl.l of the Planck reference sky, as for the CMB case. In 
this section we describe how the various components have been 
modeled. SZ sources and extra-Galactic radio sources have been 
added since Ashdown et al. (2007b). 

3.7.1. Diffuse emission 

We include synchrotron emission from free electrons spiraling 
around the Galactic magnetic field and bremsstrahlung emit- 
ted by electrons scattering onto hydrogen ions. We also include 
the emission from thermal dust grains; although subdominant 
with respect to the other components at intermediate and high 
Galactic latitudes, the brighest dusty emission regions across the 
Galactic plane are still relevant at 30 GHz. The total intensity 
information on these components is obtained from non-PLANCK 
frequencies, 408 MHz and 3000 GHz for synchrotron and dust, 



respectively, as well as Ha regions tracing the bremsstrahlung. 
No comparable all-sky information exists for the linear polar- 
ization component. The latter has been simulated by exploiting 
data at low and intermediate latitudes in the radio and microwave 
bands (see Ashdown et al. 2007b for details). 

3.7.2. Extra-galactic radio sources 

Emission from unresolved extra-Galactic radio sources has been 
obtained from existing catalogues as well as models, extrapolat- 
ing to 30 GHz. The input cata logues were the NRAO VLA Sky 
Survey (NVSS, Condon et al. |1998l l and the Sydney University 
Molonglo Sky Survey (SUMSS, Mauch et al. l2003l at 1.4 GHz 
and 0.843 GHz, respectively, which cover only part of the sky, 
as well as the Parkes-MIT-NRAO (PMN, Wright et al. |1996ab 
survey at 4.85 GHz, which covers the entire sky except for tiny 
regions around the poles. The catalogues were combined by 
degrading and smoothing the higher resolution observations to 
match those of the lower resolution surveys. To avoid double 
counting of background sources, the average flux of the NVSS 
and SUMSS surveys was evaluated after the removal of the one 
at 4.85 GHz. That average flux was then subtracted from the 
summed 4.85 GHz sources in these higher resolution surveys. 
In order to obtain a uniform map and account for the fact that 
the NVSS and SUMSS have only partial sky coverage, sources 
were copied randomly into the survey gaps from other regions 
until the mean surface density as a function of the ~ 1 GHz flux 
was equal to the overall mean down to 5 mJy. Note that the per- 
centage of simulated sources is small(~4%) and mostly located 
in the Galactic plane. 

The frequency extrapolation proceeds as follows. Sources 
were divided into two classes according to their spectral index 
a, where the flux S scales as v~ a : flat spectrum with a < 0.5, 
and steep spectrum with a > 0.5. Sources measured at a sin- 
gle frequency were assigned randomly to a class, with a drawn 
from two gaussian probability distributions, one for each class, 
with mean and variance estimated from the sample of sources 
with flux measurements at two frequencies. In the extrapolation 
at 30 GHz, corrections to the power law approximation were 



M. A. J. Ashdown et al.: Making Maps from Planck LFI 30GHz Data with Asymmetric Beams and Cooler Noise 



7 



accounted for by including the multifrequency data from the 
WMAP (Bennett et al. |2003)) l in order to derive distributions 
of differences, 6a, between spectral indices above and below 
20 GHz. For polarization, the polarization angle was drawn ran- 
domly from a fiat prior over the [0, n] interval, while the polar- 
ization percentage was drawn from a probability distribution de- 
rived from observations concerning the fiat and steep spectrum 
sources at 20 GHz. 

In the generation of TODs, sources with a flux above 
200 mJy were treated through the point source convolver code 
within the Level-S package; the remaining sources were used to 
generate a sky map that was added to the diffuse emissions. 

3.7.3. Sunyaev-Zel'dovich effect from Galaxy clusters 

We used the Monte Carlo simulation package developed by 
Melin et al. (2006) to generate the SZ cluster catalog in a 
ACDM cosmology (ACDM = cold dark matter with dark en- 
ergy). Cluster mass M and redshift z were sampled according 
to the mass function by Jenkins et al. (12001b . and we placed 
the clusters uniformly on the sky, ignoring any spatial correla- 
tions. The primordial normalization, parametrized by as, the av- 
erage mass variance within spheres of 8/z Mpc, was chosen to 
be 0.99. We normalized the temperature-mass relation follow- 
ing Pierpaoli et al. (2003 ) to match the local X-ray temperature 
function, T* = 1.3 keV, and cut the input catalog at 10 14 solar 
masses. 

The simulation assigns velocities to the cluster halos from 
the velocity distribution with variance calculated according to 
linear theory. These velocities could be used to calculate the po- 
larized SZ signal, although this feature was not implemented in 
this work. The SZ simulations in this paper are therefore unpo- 
larized. Future simulations will include SZ polarization. 

We attributed to each cluster halo an isothermal /3-model gas 
profile at the temperature given by our adopted T-M relation, 
(see Melin et al. l2006l for details). We fixed yS = 2/3 and the core 
radius of each cluster to r c - 0.1r v ; r , i.e., one tenth of the virial 
radius r v ; r ; the latter was calculated using the spherical collapse 
model. The remaining quantity is the total gas mass (or central 
density), which we determined by setting the gas mass fraction 
/gas = 0.9 • Ob/£^m (baryons and total matter). Thus the cata- 
logue is characterized by mass, redshift, position on the sky, gas 
temperature, and density profile. From this information we cal- 
culate the total integrated SZ flux density, S v , at the observation 
frequency, and then divide the catalog at lOmJy into a set of 
bright and faint sources. The bright catalog contained ~20,000 
sources that were used by the point source convolver code in the 
Level-S package to generate beam smoothed SZ point-sources 
in the TODs. We combined the catalog of fainter clusters into a 
sky map that was added to the other diffuse emissions. 

4. Mapmaking codes 

Two characteristics of mapmaking codes are important. One is 
accuracy, that is, how close a given code comes to recovering 
the input sky signal in the presence of noise and other mission 
and instrumental effects. The other is resources required, that is, 
how much processor time, input/output time, memory, and disk 
space are required to produce the map. 

Ideally, accuracy could be maximized and resource require- 
ments minimized in one and the same code. Not surprisingly, 
this is not the case. However, one can imagine different regimes 
of mapmaking, with different requirements. On the one hand, 



high-accuracy will be of paramount importance for the Planck 
legacy maps. Because such maps need be produced infrequently, 
the code can be quite demanding of resources if necessary. On 
the other hand, resources required will be critical in the interme- 
diate steps of the Planck data analysis (e.g., in systematics de- 
tection, understanding, and removal), where a great many maps 
must be made, and where Monte Carlo methods will be needed 
to characterise noise, errors, and uncertainties. 

We used mapmaking codes of two basic types, "destripers", 
and "optimal" codes (sometimes called generalized least squares 
or GLS codes, notwithstanding the fact that destriping codes also 
solve GLS equations). Key features of the mapmaking codes are 
summarized in Table 2. 

MADmap, MapCUMBA, and ROMA employ optimal algo- 
rithms, in the sense that they compute the minimum-variance 
map for Gaussian-distributed, stationary detector noise (see 
Wright [T996bl Borril l [19991 Dore et al. EO Qll Na toli et al. l20T)Tl 
Yvon & Mayet 2005 , and de Gasperis et al. 2005 for earlier work 
on optimal mapmaking). The three codes operate from similar 
principles and solve the GLS mapmaking equation efficiently 
using iterative conjugate gradient descent and fast Fourier trans- 
form (FFT) techniques. To be accurate, these codes require a 
good estimate of the power spectrum of noise fluctuations. 

Springtide and Madam employ destriping algorithms. They 
remove low-frequency correlated noise from the TOD by fit- 
ting a sequence of constant offsets or "baselines" to the 
data, subtracting the fitted offsets from the TOD, and binning 
the map from the cleaned TOD (see Burigana et al. 119971 
Delabrouille [T99"8l Maino et al. 1999, 2002, Revenu et al. l2000l 
Sbarra et al. [20031 Keihiinen et al. 2004, 2005] for earlier work 
on destriping and Efstathiou 2005, 2007 for destriping errors). 
As we will see, baseline length is a key parameter for mapmak- 
ing codes. Baseline length is adjustable in destriping codes. In 
the short baseline limit, the destriping algorithm (with priors on 
the low-frequency noise) is equivalent to the optimal algorithm. 
Similarly, optimal codes may be considered as destripers with a 
baseline given by the detector sampling rate. 

The destriper Springtide operates on scanning rings. First, 
it compresses the data by binning them in 1-hour ring maps. It 
then solves for and subtracts an offset for each ring map, and 
constructs the final output map. Due to compression of data to 
rings, Springtide can run in small memory, but its long base- 
lines (1 hour) leave larger residuals in the map at small angular 
scales. A recent feature allows Springtide to compute the hour- 
long ring maps using a one minute baseline destriper, then the 
final output is constructed as before. This double-destriping im- 
proves the maps at the cost of longer runtime. 

MADmap, Springtide, and Madam can use compressed 
pointing information, meaning that they can interpolate the de- 
tector orientation from the sparsely sampled (1 Hz) satellite atti- 
tude measurements. The other mapmaking codes require the full 
set of detector pointings sampled at the detector sampling rate. 
The main benefits of the compressed pointing are significant sav- 
ings in disk space and I/O. 

Full descriptions of the codes have been given in our previ- 
ous papers (Poutanen et al. 120061 Ashdown et al. |2007al [2007bt . 
Changes from previous versions are detailed below. 

4.1. Madam 

Madam is a destriping code with a noise filter. The user has the 
option of turning the noise filter off, in which case no prior in- 
formation on noise properties is used. Mapmaking with Madam 
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TABLE 2 

Features 



Code Madam MADmap MapCUMBA ROMA Springtide 



Algorithm Destriping Optimal Optimal Optimal Destriping 

Noise estimate needed Optional" Yes Yes Yes No 

Baselines > 1 s . . . b . . , b . . . b > 1 min 

Compressed pointing Yes Yes No No Yes 

Shifted pointing 1 -' Yes Yes Yes No Yes 

Small memory mode Yes Yes No No Yes 

Used in DPC d pipeline LFI ... HFI LFI HFI 



a Noise estimate is needed for short (< 1 min) baselines. 

b Optimal codes may be considered as destripers with a baseline given by the detector sampling rate. 

To correct the pointing shift caused by the sample integration. 

d Data Processing Center. 



for the case of noise filter turned off is discussed in Keihanen et 
al. 120081 

The baseline length is a key input parameter in Madam. The 
shorter the baseline, the more accurate are the output maps. It 
can be shown theoretically that when the baseline length ap- 
proaches the inverse of the sampling frequency, the output map 
approaches the optimal result. 

A number of improvements have been made to Madam since 
the Paris round of simulation (Ashdown et al.|2007b). The code 
constructs the detector pointing from satellite pointing, saving 
disk space and I/O time. In case two detectors have identical 
pointing, as is the case for a pair of LFI detectors sharing a horn 
antenna, pointing is stored only once, dropping the memory re- 
quirement to half. 

Further, the code uses a lossless compression algorithm 
which greatly decreases the memory consumption at long base- 
lines. Madam also allows a "split-mode", where the data are first 
destriped in small chunks (e.g., 1 month) using short baselines. 
These chunks are then combined and re-destriped, using longer 
baselines. The split-mode decreases memory consumption sub- 
stantially. The cost is that run time increases and map quality de- 
creases somewhat as compared to the standard mode. The split- 
mode can be used in many ways. One may for instance destripe 
data from 12 detectors in 3 parts, each consisting of data from 
four detectors. 

With these improvements, Madam offers wide flexibility in 
terms of computational resources used. The most accurate maps 
are obtained with a short baseline length and fitting all the data 
simultaneously. This alternative requires the maximum memory. 
The memory requirement can be reduced either by choosing a 
longer baseline, or by using the split-mode. 

The Madam maps of this study were destriped using a single 
set of short baselines (i.e., the split-mode was not used). Unless 
otherwise noted the baseline length was 1 .2 s. 

The earlier Polar code (Ashdown et al. I2007al I2007bl) . cor- 
responding to Madam with noise filter turned off, is now merged 
into Madam. 

4.2. Springtide 

A number of improvements have been made to Springtide since 
the work reported in Ashdown et al. ( |2007al l and (2007b). 

Springtide now uses the M3 data abstraction library to read 
TOD and pointing (see section |4~4l i. Instead of reading the detec- 
tor pointing information from disk, M3 can use the Generalised 
and Compressed Pointing (GCP) library to perform an on-the- 



fly calculation of the positions of the detectors from the satellite 
attitude data. 

Springtide is now capable of making maps at a number of 
resolutions in the same run. This requires that a hierarchical pix- 
elation such as HEALPix be used for the maps. The destriping 
must be performed at a sufficiently high resolution so that the 
sky signal is approximately constant across a pixel. Once the 
offsets describing the low-frequency noise are subtracted from 
the rings, they can be binned to make the output map at any res- 
olution equal to or lower than that used for the destriping. 

4.3. MapCUMBA 

The current version (2.2) has been modified in several ways rel- 
evant to this study. 

While the pointing information is generally provided in 
spherical coordinates, the mapmaking algorithm only requires 
the pixel indexing. To reduce the memory expense of storing 
both forms of pointing information while mapping one to the 
other, the pointing is read from disk into a small buffer (whose 
size can be adjusted by the user) and is immediately mapped into 
the final pixel index stream. The drawback of this scheme is that 
it is generally much faster to read the same amount of data from 
disk in one piece than in several small pieces. In the configura- 
tion tested, the I/O buffer has to be kept larger than 10 7 samples 
for the I/O not to dominate the total run time. 

Since the preconditioned conjugate gradient (PCG) al- 
gorithm used in MapCUMBA involves repeated overlap-add 
Fourier transforms of fixed length, it is beneficial to use highly 
optimized FFT algorithm such as fftw-fl We found the ID 
Fourier transform offered by fftw-3.0.0, with a 'measured' plan 
selection algorithm and a length of 262144, to be twice as fast as 
the one implemented in fftw-2. 1 .5, offsetting both the more cum- 
bersome interface of fftw-3, and the overhead associated with the 
'measured' plan selection over the 'estimated' one. 

4.4. MADmap 

MADmap is the optimal mapmaking component of the 
MADCAP suite of tools, specifically designed to analyse large 
CMB data sets on the most massively parallel high performance 
computers. Recent refinements to MADmap include options to 
reduce the memory requirement (at the cost of some additional 
computations), and to improve the computational efficiency by 
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choosing the distribution of the time ordered data over the pro- 
cessors to match the requirements of a particular analysis. Like 
Springtide, MADmap uses the M3 data abstraction and the GCP 
libraries. 

M3 allows an applications programmer to make a request to 
read a data subset that is independent of the file format of the 
data and the way the data are distributed across files. In addi- 
tion it supports "virtual files", which do not exist on disk and 
whose data are constructed on the fly — specifically used here 
by MADmap to construct the inverse time-time noise correla- 
tion functions from a spectral parametrization of the noise. This 
abstraction is mediated through the use of an XML description 
of the data called a run configuration file (or runConfig), which 
also provides a convenient way of ensuring that exactly the same 
analysis is executed by different applications — MADmap and 
Springtide here. 

GCP provides a way to reduce the disk space and 10 require- 
ments of mapmaking codes. Instead of storing the explicit point- 
ing solution for every sample of every detector, we store only the 
pointing of the satellite (generalized) every second (compressed) 
and reconstruct via M3 the full pointing for a particular set of 
samples for a particular detector through on-the-fly interpola- 
tion and translation only when it is requested by the application. 
In this analysis — mapping only the four slowest-sampled of 72 
Planck detectors, comprising a little over 1% of the data — the 
use of the GCP library reduced the disk space and IO require- 
ments for MADmap and Springtide from 92 to 1.5 GB. 

MADmap uses the PCG algorithm to solve the GLS prob- 
lem. This requires applying the pointing matrix and its transpose 
once on each iteration of the PCG routine. When running in high 
memory mode the sparse satellite pointing is expanded once us- 
ing GCP, and the entire portion of the pointing matrix for the 
time samples assigned to a processor is stored in its memory (in 
packed sparse form) to be reused in each PCG iteration. When 
running in low memory mode only small portions of the point- 
ing matrix that fit into a buffer are expanded in sequence and the 
buffer sized pointing matrix is used and then overwritten. This 
exchanges memory used for computation time spent within the 
GCP library, but this allows for the analysis of very large data 
sets on systems where the memory requirements would other- 
wise be prohibitive. The GCP library has been optimized to be 
very computationally efficient, including the use of vectorized 
math libraries, and typically consumes about 1/3 of the run time 
in a MADmap job in low memory mode. In low memory mode 
the memory consumption scales with the number of pixels, es- 
sentially independent of the (very much larger) number of time 
samples. 

The other feature recently added to MADmap is an alter- 
native distribution of time ordered data over the processors. The 
only distribution available in previous versions of MADmap was 
to concatenate all the detectors' time streams into a single vec- 
tor and distribute this over the processors so that each processor 
analyses the same number of contiguous detector samples. This 
distribution is still an option in MADmap, but there is now an 
alternative which saves memory and cycles in certain circum- 
stances by reducing the amount of compressed pointing data that 
each processor must calculate and store. In the alternative dis- 
tribution each processor analyses data for a distinct interval of 
time for all of the detector samples that occur in that time inter- 
val. These time intervals are chosen so that each processor has 
the same total number of samples (regardless of gaps in detector 
data). This implies that each processor analyses data from ev- 
ery detector. Note that each processor stores a distinct portion of 
the compressed pointing (modulo small overlaps to account for 
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noise correlations). This distribution will also allow future ver- 
sions of MADmap to include the analysis of inter-channel noise 
correlations. The runs described in this paper were all done with 
the original concatenated distribution for time ordered data. 

4.5. ROMA 

ROMA is now stable at version 5.1 (the same as was employed 
for Ashdown et al. d2007al) . d2007bl) ). which makes use of fftw- 
3.1.2. Nonetheless, a few minor improvements have been per- 
formed to optimize speed (by tuning some fftw parameters) and 
memory usage. In addition, a new I/O module has been devel- 
oped for the sake of the present simulations that allows time- 
lines containing different sky components to be read in quickly 
and then mixed together. 

5. Results of mapmaking 

In this section we quantify the results of our mapmaking exer- 
cise. Our main goal was to examine the effects of detector main 
beams. To compare the maps and the effects of systematics upon 
them, we define three auxilliary maps: 

- input map represents the true sky. In some cases it contains 
the CMB alone; other cases it includes the dipole and fore- 
grounds as well. The CMB part of the input map contains no 
B-mode power. 

- smoothed input map is the input map smoothed with an ax- 
ially symmetric Gaussian beam and a pixel window func- 
tion^ 

- binned noiseless map refers to the (/, Q, U) map obtained 
by summing up the noiseless time-ordered data, accounting 
for the detector orientation. This is the best map of an ideal 
noiseless situation. 

The binned noiseless map (m B ) is produced from the noise- 
less TOD (s) as 

m B = (P r P) ' P r s, (2) 

where P is the pointing matrix, which describes the linear com- 
bination coefficients for the (/, Q, U) pixel triplet to produce a 
sample of the observed TOD. Each row of the pointing matrix 
has three non-zero elements. 

The temperature of a pixel of the smoothed input map is an 
integral of the (beam smoothed) sky temperature over the pixel 
area. In the binned noiseless map the corresponding pixel tem- 
perature is not a perfect integral, but a mean of the observations 
falling in that pixel. This pixel sampling is not as uniform as the 
integration in the smoothed input map. Therefore the difference 
of the binned noiseless map and the smoothed input map has 
some pixel scale power due to this difference in pixel sampling. 
We call this difference pixelization error in this paper. For asym- 
metric beams there is the additional effect that different observa- 
tions centered on the same pixel may fall on it with a different 
orientation of the beam, resulting in a different measured signal. 

For CMB, dipole, and foreground emissions we made four 
different simulated TODs, depending on whether the beams were 
axially symmetric or asymmetric Gaussians, and whether the 
sample integration was on or off (see Sect. [21). 

6 For the axially symmetric beam, we used FWHM = 32'. 1865, which 
is the width of the symmetric beams of this study (see Sect. l3.2b . For the 
pixel window we used the HEALPix pixel window function of N S id e = 
512 pixel size (Gorski et al. 2005bl. 
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The widths of the symmetric beams were identical in all four 
LFI 30 GHz detectors, whereas the widths and orientations of 
the asymmetric beams were different. The difference of the de- 
tector beam responses is called the beam mismatch. For the case 
of sample integration, the time stamps of the detector pointings 
were assigned to the middle of the sample integration interval. 
When the sample integration was turned off, the observations of 
the sky signal were considered instantaneous and the timings of 
the detector pointings coincided with them. 

The mapmaking methods discussed here utilize the detector 
pointing information only to the accuracy given by the output 
map pixel size. The methods use the pointing matrix P (see Eq. 
(f2]i) to encode the pointings of the detector beam centers and 
the directions of their polarization sensitive axes. None of our 
mapmaking codes makes an attempt to remove the beam convo- 
lution from its output map. Therefore an output map pixel is con- 
volved with its own specific response {effective beam), i.e., the 
mean of the beams (accounting for their orientations) falling in 
that pixel. Recently, deconvolution mapmaking algorithms have 
been developed that can produce maps in which the smoothing 
of the beam response has been deconvolved. These methods lead 
to maps that approximate the true sky (Burigana & Saez 2003 
Armitage & Wandelt 120041 Harrison et al. 120081 ). These meth- 
ods are sufficiently different from the ones considered here that 
different methods of comparison must be used, and we therefore 
do not include detailed description of them in this paper. 

An output map of a mapmaking code can be considered as 
a sum of three components: the binned noiseless map, the resid- 
ual noise map; and an error map that arises from the small- 
scale (subpixel) signal structure that couples to the output map 
through the mapmaking (Poutanen et al. 2006 Ashdown et 
al. 1200 7b). We call the last error map the signal error map. We 
use the term residual map, when we refer to the sum of the sig- 
nal error and residual noise maps. Because the signal error arises 
from the signal gradients inside the output map pixels, smaller 
pixel size (i.e., higher map resolution) leads to a smaller signal 
error. Our earlier studies have shown that for a typical Planck 
map (e.g., A^ide = 512 or smaller pixels), the signal error is a 
tiny effect compared to the CMB signal itself or to the residual 
noise (Poutanen et al. 120061 Ashdown et al. |2007bl > . Therefore 
a Planck signal map is nearly the same as the corresponding 
binned noiseless map. It is a common characteristic of all our 
maps. 

Generally, the binned noiseless map contribution in the out- 
put map is (cf. Eq. (f2]l) 

m B = (P^C^P) - ' P r C;'s, (3) 

where C„ is the diagonal time-domain covariance matrix of the 
detector white noise floor. Its diagonal elements will not equal in 
general; for example, the white noise RMS of the detectors can 
be different. However, in this study all detectors are assumed 
to have the same white noise and therefore the matrix C„ can 
be ignored and Eq. (f2]i gives the correct binned noiseless map 
contribution. 

In this study we assume that the residual noise map contains 
residues of the uncorrelated (white) and correlated (1//) instru- 
ment noise and the residues of the sorption cooler fluctuations. 

The (/, Q, U) maps we made in this study were pixelized at 
A^side = 512. At this resolution every pixel was observed (full sky 
maps) and their polarization directions were well sampled. The 
rcond's of the 3x3 N b s matrices were larger than 0.3q 

7 The quantity rcond, the reciprocal of the condition number, is the 
ratio of the absolute values of the smallest and largest eigenvalue of the 



Unless otherwise noted the maps are presented in ecliptic co- 
ordinates and thermodynamic (CMB) microkelvins. The units of 
angular power spectra are thermodynamic microkelvins squared. 

To demonstrate the effects of beams in our maps, we made 
maps from the noiseless TODs containing CMB, dipole, and 
foreground emissions. We show the Madam maps as an example 
in Fig. [5] Corresponding maps of the other mapmaking codes 
would look similar. We can hardly see any differences between 
the noiseless output map and the smoothed input map (see the 
two upper rows of Fig. 01. To reveal the differences, we sub- 
tract the smoothed input map from the output map (bottom two 
rows of Fig. 0. The beam window functions of the symmetric 
and asymmetric beams differ mainly at high-^ (see Sect. I5.1.2I I. 
which makes mainly small angular (pixel) scale differences in 
the maps. Therefore the difference map of the asymmetric beams 
contains more small-scale residuals than the difference map of 
the symmetric beams. The ecliptic pole regions of the sky are 
scanned in several directions, which makes the effective beam of 
the asymmetric case more symmetric and therefore closer to the 
effective beam of the symmetric case. Therefore the difference 
of the noiseless output map of the asymmetric beams and the 
smoothed input map becomes small in the vicinity of the eclip- 
tic poles (see the light green areas of the difference map of the 
third row of Fig. |5). The angular diameter of these areas is -30°. 

Some point source residues are visible in the temperature dif- 
ference maps of Fig.|5](see especially the lower left corner map). 
A likely reason for these residues is the difference in how the 
pixel areas are sampled in the output maps and in the smoothed 
input map. In the latter map the pixel temperature is an integral 
of the sky temperature over the pixel area, whereas in the out- 
put map the pixel temperature is an average of the observations 
falling in the pixel. The observations do not necessarily sample 
the pixel area uniformly, which leads to a different pixel temper- 
ature than in the uniform integration. 

The difference maps of Fig. [5] have some stripes that align 
with the scan paths between the ecliptic poles. These stripes are 
most noticeable close to the galactic regions in the symmetric 
case difference map. In the asymmetric case the stripes do not 
stand out from the larger pixel scale residuals. They arise from 
the fact that signal differences (gradients) inside a map pixel 
create non-zero baselines in Madam that show up as stripes in 
the Madam map. Because the signal gradients are largest in the 
galactic regions, the stripes are strongest there. All our map- 
making codes produce such signal errors, which are stronger in 
the optimal codes than in the destripers (Poutanen et al. 2006 
Ashdown et al.l2007bj. 

Another way to see the effects of the beams in the maps is to 
examine how point sources show up in the maps. The image of 
a point source in the map shows the effective beam at that loca- 
tion of the sky. Figs.|6]and|7]show two such point sources of the 
noiseless Madam temperature maps. Fig. [6] is a patch from the 
vicinity of the ecliptic plane. There the scanning is mainly in one 
direction and therefore the difference in ellipticities of the effec- 
tive beams (of the symmetric and asymmetric cases) is clearly 
visible. Fig.|7]shows a similar comparison near the south eclip- 
tic pole. Here the wide range of the scanning directions makes 
the effective beams more symmetric. The effective beams of the 
symmetric and asymmetric cases are now more alike, but we can 
still detect some difference in their ellipticities. 



3x3 N (,s matrix of a pixel. The matrix P r P is block-diagonal, made 
up of these N bs matrices. For a set of polarized detectors with identical 
noise spectra (like the LFI 30 GHz detectors of this study) rcond is < 0.5 
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Fig. 5. Madam maps made from simulated noiseless TODs. Maps on the left are temperature maps. Those on the right give the mag- 
nitude of the polarization vector (P = V2 2 + U 1 ). All maps contain CMB, dipole and foreground emissions. Top row: Smoothed 
input maps. Second row: Madam maps. These maps include the effects of asymmetric beams and sample integration. Third row: 
Difference of the above Madam map and the smoothed input map. Bottom row: Same as above, but for symmetric beams. The maps 
of the top and second rows are in ecliptic coordinates, whereas the difference maps are in galactic coordinates. The latter coordinates 
were chosen to give a clearer view to the ecliptic pole regions (see the light green areas of the asymmetric difference map). Small 
patches of noiseless Madam temperature maps are shown in Figs. [6] and [7] 



In our simulations beams and sample integration distort the maps, independent of any particular mapmaking algorithm. We 
TOD before the instrument noise is added. Independent of the examine the effects of the beams and sample integration on the 
noise, these effects are best explored in the binned noiseless binned noiseless maps in Sect. 15. II 
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Fig. 6. Effect of beams in the point source observations. These 
plots show 3° x 3° patches of the noiseless Madam temperature 
maps (from Fig. ID, zoomed to the vicinity of the ecliptic plane. 
The location of the patch is given as a square box in the top left 
map of Fig. [5] 
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Fig. 7. Similar plots as in Fig. [6] but now we have zoomed to a 
point source near the south ecliptic pole (in (A,fi) = (255°, -81°) 
of ecliptic longitude and latitude). The size of the patch is 3°x3°. 



Because the beams and sample integration affect the signal 
gradients of the observations, they have an impact on the signal 
error too. We examine these effects in Sect. 15. 21 and compare dif- 
ferences between mapmaking codes. It is only through the signal 
error that the effects of beams and sample integration show up 
differently in the maps. The binned noiseless map contribution 
that is also affected by beams and sampling stays the same in 
all maps. In our simulations beams should have no effect in the 
residual noise maps. Because there is a half a sample timing off- 
set between the detector pointings of sampling on and off cases, 
the hit count maps of these two cases differ slightly, which shows 
up as small differences in the residual noise maps. 

Finally we discuss sorption cooler fluctuations and detector 
pointing errors and assess their impacts in the maps. 

5.1. Binned noiseless maps 

The effects of beams and sample integration were examined 
from the binned (/, Q, U) maps made from the four noiseless 
TODs of different beam and sampling cases. The TODs con- 
tained only CMB. The angular power spectra of the four maps 



are shown in Fig. [8] The TT angular power spectra of the sym- 
metric and asymmetric beams differ mainly at large multipoles 
({ > 600). This is a result of the different effective beam window 
functions of these two cases. 

The EE spectra of symmetric and asymmetric beams are dif- 
ferent at both intermediate (I > 100) and high multipoles. Due 
to the beam mismatch the EE spectrum of the asymmetric beam 
case is influenced by the cross-coupling from the temperature. 
A similar cross-coupling does not occur in the symmetric beam 
case (no beam mismatch). Therefore the EE spectra behave dif- 
ferently than the TT spectra. We discuss these issues more in 
Sects. IBTTTl and I5.1.2l and give there the explanation of the be- 
haviors of the TT and EE spectra. 

Sample integration is effectively a low-pass filter in the TOD 
domain. Therefore it introduces an extra spectral smoothing that 
removes some small-scale signal power from the map. This ef- 
fect is just barely visible in Fig. [8] where the TT and EE spec- 
tra with "sampling" are slightly suppressed as compared to their 
"no sampling" counterparts. We discuss these effects more in 
Sect. 15X41 

The TT and EE spectra of the binned maps become flat at 
high £ (> 1000). This is a result of spectral aliasing (£ mode 
coupling) that arises from the non-uniform sampling of the pixel 
areas. The aliasing couples power from lovi-i to high-^. We see 
this effect in the maps too, where we called it the pixelization 
error. 

Because the CMB input sky contained no B-mode power, the 
signals in the BB spectra of the binned maps (see Fig. [8]l must 
arise from temperature and /s-mode polarization signals leaking 
to the B-mode. The magnitudes of the symmetric beam BB spec- 
tra are comparable to the magnitude of the high-f flat part of the 
corresponding EE spectra. This observation and the fact that the 
symmetric beams were all identical (no beam mismatch) suggest 
that the BB spectrum of the symmetric beam case mainly arises 
from the spectral aliasing from the E mode (due to pixelization 
error). The BB spectrum of the asymmetric beam case shows a 
distinct signal at multipoles between 100 and 800. We expect 
this signal to originate from temperature and E-mode polariza- 
tion signals that cross-couple to BB due to beam mismatch. We 
discuss these cross-couplings more in Sect. 15.1.11 

5.1 .1 . Beam mismatch, cross-couplings of /, Q, and U 

We can quantify the effects of beam mismatch on the maps by 
calculating power spectra. We consider the asymmetric beam 
case with no sampling A number of authors have worked on 
beam mismatch systematics and their impacts on maps and an- 
gular power spectra (Hu et al. l2003l Rosset et al. I2007I O'Dea et 
al. |20U71 Shimon et al. |2008l . 

The LFI main beams vary in width and orientation from de- 
tector to detector. These variations are fully represented in our 
asymmetric beams. (In contrast, the widths of our symmetric 
beams were the same in all detectors.) Due to beam mismatch, 
the detectors of a horn see different Stokes I; this difference ap- 
pears as an artifact in the polarization map. This is a potentially 
serious issue for a CMB experiment such as Planck, because the 
fraction of the strong temperature signal that pollutes the polar- 
ization map (denoted T — > P) may be fairly large compared to 
the weak CMB polarization signal itself. 



8 Note that for this analysis, the smoothing effect from detector sam- 
ple integration could be subsumed into the asymmetric beams, but we 
decided to ignore it for simplicity. 
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Fig. 8. Overview of the effects of beams and sample integration on the angular power spectra of the CMB maps. We computed the 
spectra of the noiseless CMB maps that were binned from the four simulated LFI 30 GHz TODs. The TODs represent different 
beam and sample integration cases. The CMB input sky contains no B mode power. Note that the plots of the cross-spectra show 
their absolute values. In this paper we define Cf Y = 1/ m =-( a J ( a J )*■ Here X,Y — T,E,B and Z* is the complex conjugate of 
Z. 



Strong foreground emission dominates the T — > P signal 
in a full sky 30 GHz map. This map is useful for CMB power 
spectrum estimation only if the pixels with strong foreground 
contribution are removed from the map (galactic cut); however, 
the T — > P signal will depend on the mask used in the cut. To 
avoid this, and to study the beam mismatch effects in a more 
general case for which our results would not depend much on 
the details of the data processing (e.g., masks), we decided to 



continue to work with the full sky noiseless maps that we binned 
from the CMB TOD. 

To give low-^ details we replotted the EE and TE spectra 
(from Fig. [8) of the binned noiseless CMB map (asymmetric 
beams & no sampling) and the smoothed input map. The replot- 
ted spectra are shown in Fig. [9] The T — > P cross-coupling is a 
significant systematic for the power spectrum measurement, as 
Fig.|9]shows. The acoustic peaks and valleys of the binned map 
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are shifted towards higher multipoles. In the EE spectrum the ef- 
fect is smaller than the noise, and may not be detectable. For the 
TE spectrum the error of the E mode polarization gets amplified 
by the large temperature signal. This results in a larger artifact 
that is visible even in the noisy TE spectrum (see Fig.[20]and the 
discussion in Sect. I5.2t . 

We designed a simple analytic model that gives a reasonably 
good description of the effects of the beam mismatch in the an- 
gular power spectra of the (7, Q, U) signal maps (Appendix [A}. 
Our model shows that the mismatch of the beams of a horn is 
relevant for T — > P. The beams can be different from horn to 
horn, but only a weak T — > P may arise if the two beams of a 
horn are identical. 

To isolate the effect of T — > P, we generated test TODs for 
all four 30 GHz radiometers. The TODs contained Stokes / of 
the CMB only. We binned a new (/, Q, U) map from these TODs. 
The polarization part of the map should contain the temperature 
leakage only. We therefore call this map the leakage map. The 
power spectra of the leakage map and the original binned map 
with the leakage map subtracted are shown in Figs. [TOl and [TTI 
Removing T —> P restores the acoustic peaks and valleys of 
the EE spectrum in their original positions. The mismatch of the 
beams is important at or below the beam scale. Therefore the 
magnitude of the T — > P leakage, and the corresponding error in 
the EE spectrum, are small at low (. 

Beam mismatches cause cross-coupling in the opposite di- 
rection too, namely E — > T and B — > T; however, because the 
E and B mode powers are small compared to the power of the 
temperature signal, these couplings have an insignificant effect 
in the maps and angular power spectra. 

The B-mode spectra (left-hand panel of Fig. [T2l show two 
notable effects. First, the cross-coupling from the temperature to 
the E- and B-mode spectra are not the same. Though similar in 
shape, the T — > E power is 2-3 orders of magnitude larger than 
the T — > B power (at I ~ 200). This discrepancy is rooted in 
the beam widths and orientations coupled with the scan pattern. 
Our analytic model (of Appendix [A} is able to predict this be- 
havior (see the right-hand panel of Fig. 112). The analytic model 
also shows that if the detector beams were axially symmetric but 
different in widths, E and B would show the same amount of 
leakage from temperature. 

Second, the B-mode spectrum shows more power than can 
be explained by T — > P alone (see the green curve of Fig. [12). 
Because the CMB sky did not contain B-mode polarization, the 
B-mode power remaining after the removal of the T — > P must 
result from the E — > B cross-coupling. The source of the dis- 
tinct E — > B signal at i < 800 is the spin-flip coupling. In 
Appendix [A] we show that the relevant quantity in E — > B is 
the sum of the beam responses of the pair of detectors sharing 
a horn. If there is a mismatch of these sums (mismatch between 
the horns), the two polarization fields of opposite spins (Q + ill 
and Q — iU) get mixed. This is the spin-flip coupling (Hu et 
al. 120031 1. Appendix lAl shows how the spin-flip coupling arises 
from the beam mismatch and how it creates E — > B. In these 
simulations the source of the mismatch of the sum responses is 
not the widths of the beams, but their orientations (see Sect. 13.2) . 
The widths of the pair of beams of a horn are different, but the 
pairs of beams have these same values in both horns. The orien- 
tations of these pairs are, however, different in the two horns. 

In spite of the fact that the CMB /s-mode polarization sig- 
nal is significantly weaker than the temperature signal, the mag- 
nitude of the E — > B signal is larger than the magnitude of 
T — > B signal. It seems that the difference of the orientations 
of the pairs of beams produces an E — » B signal that is stronger 



than the T — > B signal produced by the mismatch of widths of 
the two beams of a horn. The power transfer between the po- 
larization modes operates equally in both directions. B — > E 
cross-coupling occurs too (same coupling transfer function as in 
E — > B), but it has no effect in E, because the B-mode power 
is zero. At large multipoles (flat part of the BB spectrum at 
{ > 800) the main source of E — > B is the I mode coupling 
of the non-uniform sampling of the pixels (pixelization error). 

Finally, the T — > B and E — > B signals could be compared 
to the magnitude of the CMB B-mode signal that we might ex- 
pect to detect with LFI 30 GHz detectors. In the right-hand panel 
of Fig.Q~2]the leakage spectra are compared to a theoretical BB 
spectrum of the CMB (including lensing from E mode and corre- 
sponding to a 10% tensor-to-scalar ratio). It can be seen that the 
cross-coupling to B is small compared to this signal for £ < 300. 

5.1 .2. Effective window functions 

The effects of instrument response and data processing on the 
map can be described in terms of an effective window function, 
which we compute as the ratio of the map power spectrum to 
the input spectrum. In our simulations, beams, sample integra- 
tion, and sampling of the pixel area (pixel window function) are 
the main contributors to the effective window function. In some 
cases, to reveal more details, we compute the effective window 
functions relative to the smoothed input spectrum. 

Figure[T3]shows the TT, EE, and TE effective window func- 
tions of the binned noiseless CMB maps that we introduced ear- 
lier in this section. The Gaussian window function approxima- 
tion breaks down at high Z and the function becomes flat or starts 
to increase (as in the TT window functions). This is an effect 
of the pixelization error in the binned noiseless maps. Polarity 
changes of the TE spectra cause the sudden jumps seen in the 
TE window functions. The EE and TE window functions of 
the asymmetric beams show another non-regular characteristic. 
They deviate significantly from the regular Gaussian response. 
The source of this effect is the T — > P cross-coupling. The mag- 
nitudes of the window functions with "sampling on" are system- 
atically smaller than their "sampling off" counterparts. This is 
more clearly visible in the window functions of the symmetric 
beams. This difference in window functions comes from the ex- 
tra spectral smoothing of the sample integration. 

To examine the window functions without the effects of the 
T — > P cross-coupling seen above, we subtracted the T — > P 
leakage map from the polarization part of the binned noiseless 
CMB map and recomputed the window functions, as shown in 
Fig. [14] We did this for the asymmetric beams/no sampling case, 
because we had the leakage map for this case only. The win- 
dow functions have now more regular shapes. The TT and EE 
window functions are not identical; the TT function is slightly 
steeper than the EE function. In the case of perfectly matched 
beams the TT and EE window functions would be identicafl. 

The subtraction of the T — > P leakage map from the binned 
map does not influence the E — > B cross-coupling. In addi- 
tion to creating a B-mode polarization signal from the E'-mode 
signal, this coupling influences the original E-modt signal (see 
Sect. lA.ll of AppendixlAb. so that the TT and EE window func- 
tions become different. Fig. [T41 demonstrates that our analytical 
model is able to explain these window functions. 

The effective window function of the symmetric beams with 
no sampling is simple to compute because all detectors have the 

9 The off-diagonals of the Mueller matrix would be zero and all di- 
agonals would be identical (see Eq. \A.7\ ). 
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Fig. 9. Effects of T — > P cross-coupling arising from beam mismatches. Left panel: Replot of the EE spectra from Fig. [8] showing 
spectra of a binned noiseless CMB map (blue curve) and smoothed input map (red curve). The expected EE spectrum of the white 
noise map is indicated by the horizontal dashed line. Right panel: TE spectra of the same maps. We show the absolute values of the 
spectra. The TE spectrum does not have a noise bias, because white noise is uncorrected in the temperature and polarization maps. 
Therefore the dashed line of noise bias is missing from this plot. 
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Fig. 10. Demonstration of how removing the T — > P cross-coupling restores the acoustic peaks and valleys of the CMB in their 
original positions. Left hand panel: We show here a number of EE spectra. The gray curve shows the EE spectrum of the leakage 
map. The red and blue curves are as in Fig. [9] The green curve shows the spectrum of the binned noiseless map, where the leakage 
map has first been subtracted. Note that the green curve is nearly on top of the red curve at low multipoles (£ < 400). The red and 
green curves differ at t > 450 mainly due to differences in the beams: the red spectrum is smoothed with symmetric beams while 
the green spectrum is smoothed with asymmetric beams. At I > 1050 the green and blue spectra become flat due to spectral aliasing 
arising from pixelisation errors. T — > P becomes flat too for the same reason. The expected spectrum of the white noise map is 
indicated by the horizontal dashed line. Right hand panel: Same as left hand panel but zoomed to lower I. The red and green curves 
are nearly on top of each other. 



same beamwidths. If we assume that the scanning is simply from 
pole to pole and all detector beams have fixed orientations rela- 
tive to the local meridian, the effective TT window function of 
the asymmetric beams (no sampling) can be estimated in the fol- 
lowing way: 



3. Compute the mean b( m over the radiometer beams: 

bin, = jfc ■ l£j b L- FOT LFI 30 GHZ N ** = 4. 

4. Compute the effective beam window function as 

^1 ~ 2f+l ' 2m \blm\ ■ 



Denote by Bk(6, <p) the beam response of the radiometer k of 
the LFI focalplane (see Fig.|2]i. The beams are normalized to 

j 4n B k (o,<p)dn = i. 

Compute the coefficients (b k £m ) of the spherical harmonic ex- 
pansion of the beam B^O, <p). 



If this window function is convolved with the HEALPix 
pixel window function, the result explains well the TT window 
function of Fig. [14] 
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Fig. 11. Restoration of the acoustic peaks and valleys in the TE spectrum. Left hand panel: Similar to the left hand panel of Fig. [TO] 
but for TE spectra. Note that we plot the absolute value of the spectra. Right hand side: Same as the left hand panel, but we plot 
now 1(1 + \)C T e E I2n and zoom to lower multipoles. The red and green curves are nearly on top of each other. 
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Fig. 12. Left hand panel: Effects of T — > P and E — > B cross-couplings on the BB spectra. The blue curve is the spectrum of 
the binned noiseless CMB map. The gray curve is the BB spectrum of the T — > P leakage map. Because the input CMB sky had 
zero B-mode power, the remaining BB power after the subtraction of the leakage map is mostly cross-coupling from the E'-mode 
polarization. This is shown by the green E — > B curve. The expected spectrum of a white noise map is indicated by the horizontal 
dashed line. For comparison, we show a theoretical B-mode spectrum corresponding to a 10% tensor-to-scalar ratio and including 
lensing from E. Right hand panel: EE and BB spectra of the T — > P leakage map (red and blue curves, respectively) and the BB 
spectrum of the E — > B coupling (green curve). They are the same as the gray and green curves of the left-hand panel of this figure 
and the gray curve of the left-hand panel of Fig. [10] The thin black lines are the predictions of these signals computed from the 
analytical model (Appendix [A). Flattening at high £ results from non-uniform pixel sampling. This effect is not included in our 
model, however, so the model spectra do not become flat. 



5.1.3. Correction of beam mismatch effects 



For temperature observations, the aj m of the sky get convolved 
with a beam, which is fully described by one complex num- 
ber for every I and m. For polarization, where we need three 
complex quantities to describe the sky signal (aj^' 8 ), the beam 
in general is a complex 3x3 matrix (for every { and m). Its 
non-diagonal elements, which arise from the beam mismatch, 
are responsible for the cross-couplings of temperature and po- 
larization. The leakage map approach that we used earlier to re- 
move the T — > P effects is not applicable in real experiments 
where Stokes /-only timelines cannot be constructed indepen- 



dently. For real experiments, more practical methods to correct 
beam effects are required. 

Mapmaking methods that address beam convolution prop- 
erly have been proposed. Deconvolution mapmaking with a 
proper treatment of the detector beams is a method to produce 
maps free from cross-couplings of T, E, and B. Two imple- 
mentions of this method have been introduced (Armitage and 
Wandelt I2004I Harrison et al. I2008I I. Both have shown results 
indicating that they may be computationally practical for the 
lower-resolution Planck channels (in the LFI). Another map- 
domain method able to correct for beam effects is the FICSBell 
approach (Hivon et. al. 2008), in which asymmetries of the main 
beam are treated as small perturbations from an axially symmet- 
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Fig. 13. TT, EE, and TE effective window functions, computed as ratios of the angular power spectra of the binned noiseless CMB 
maps and the input CMB. Because the Z?-mode power of the CMB was zero, no BB window function could be computed. For 
comparison we show the total response of the axially symmetric Gaussian beam (dashed curve; FWHM = 32'1865 and A^ide =512 
HEALPix pixel window function). This is the window function that we applied to the input spectrum to obtain the smoothed input 
spectrum. 



ric Gaussian beam, and are averaged over each pixel taking into 
account the orientation of the detector beams at each visit of that 
pixel. 

In Appendix [A] we developed an analytical model to predict 
the effects of the beams in the angular power spectra of the CMB 
maps. This model was inverted (also in Appendix [At to turn 
it to a correction method. It can deconvolve the effects of the 
asymmetric beams from the angular power spectrum and return 
a spectrum that is an approximation of the spectrum of the input 
sky. Our method is based on a number of simplifying assump- 
tions that limit its accuracy in real experiments; however, we can 
use it to compute coarse corrections that can be improved with, 
e.g., Monte Carlo simulations. The correction capability of our 
method is demonstrated in Fig. [20l 

5.1.4. Sample integration effects 

We recomputed the effective TT window functions of Fig. [13] 
but this time we divided the angular power spectra of the binned 
noiseless CMB maps with the smoothed input spectrum (instead 
of the input spectrum). The recomputed effective TT window 
functions are shown in Fig. [15] A pairwise comparison of the 
"sampling" and "no sampling" window functions (of the same 
beam type) reveals the effect of the sample integration in the an- 



gular power spectrum of a CMB map. The extra spectral smooth- 
ing due to the sample integration can be clearly seen. 

Let us consider an axially symmetric Gaussian beam and its 
elongation in the direction of the scan. The window function of 
the initial beam (before elongation) is 



Bi = e 



(4) 



where cr = FWHM/ V8hi20 This window function operates 
in the angular power spectrum domain. The beam <x stays in 
its original value (cr = cr ) in the perpendicular direction of 
the scan. Along the scan the cr gets modified to (Burigana et 
al. [200Tb 



'eff 



12 



(5) 



Here A6* s = 27r/ sp i n // samp i e , the angle through which the beam 
center pointing rotates during a detector sample time. For the 
nominal satellite spin rate (/ sp i n = 1 rpm) the ellipticity (cr e ff /<xo) 
of the elongated beam is 1.027. The ellipticity produced by 
the scanning is significantly smaller than the ellipticities of our 
asymmetric beams (-1.35). 



10 In this study the FWHM of the axially symmetric LFI 30 GHz 
beams is 32' 1865. 
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Fig. 14. TT, EE, and TE effective window functions after the T — > P leakage map has been removed from the binned noiseless CMB 
map (asymmetric beam & no sampling). Top row: The TT window function is the ratio of the power spectra of the binned noiseless 
CMB map and the CMB input map. The EE and TE window functions were calculated similarly, except that Q and U leakage maps 
were first subtracted from the original binned noiseless map. Symmetric Gaussian responses were separately fitted to the TT and EE 
window functions (in range ( = . . .500). The fit FWHM were 32.32 (TT) and 31. '94 (EE). Bottom panel: The TT and EE window 
functions of the top left panel were divided by the dashed theoretical window function. The resulting ratios are shown in this plot. 
Thin black curves are window functions that we computed from our analytical model for this case (see Appendix IA1. 



The geometric mean of the cr's of the elongated beam is 



cr e = cr cr eff = cr cr () + 



2 



A6i 



(6) 



12 " u ' 24 

The last form is an approximation that could be made because 

^j- <K o-y. Eq. (O suggests that the effect of the sample integra- 
tion in the angular power spectra of the maps could be approxi- 
mated by a symmetric Gaussian window function 

where cr s = A8 S / a/24. In our simulations this cr s corresponds to 
FWHM = 5'32. We compare the predictions of this model to the 
actual window functions of our simulation in Fig. [15] The com- 
parison shows that the accuracy of our simple model is good in 
the symmetric beam case but somewhat worse for the asymmet- 
ric beam case. 

5.2. Residual maps 



binned noiseless map. Smaller residuals imply smaller mapmak- 
ing errors. We examined the RMS of the residuals, and computed 
power spectra to study scale dependence. We are interested in 
anisotropy; the mean sky temperature is irrelevant. Therefore, 
whenever we calculated a map RMS, we subtracted the mean of 
the observed pixels from the map before squaring. The RMS of 
a map was always calculated over the observed pixels. 

Fig. [16] shows the RMS of the residual temperature maps 
for some of our mapmaking codes and for a number of de- 
striper (Madam) baseline lengths. The data of Fig. [16] were 
derived from HEALPix N s id e =512 maps. At that pixel size 
the instrument noise is the dominant contributor in the resid- 
ual mapf0(Poutanen et al. [20061 Ashdown et al. I2007bl >. The 
optimal codes consistently deliver output maps with the small- 
est residuals, which were nearly the same between the several 
codes. Madam with short uniform baselines (e.g., 1.2 s) produces 
maps with residuals essentially the same as those of the optimal 
codes, at considerably lower computational cost. Destriper maps 
with longer baselines showed larger residuals. The RMS of the 



To compare mapmaking codes, we constructed residual maps, 
specifically the difference between the output map and the 



11 Mapmaking errors arising from the subpixel signal structure typi- 
cally show up in maps with considerably larger pixel size. 



M. A. J. Ashdown et al.: Making Maps from Planck LFI 30GHz Data with Asymmetric Beams and Cooler Noise 



19 




10 1 10 2 10 3 



Multipole I 

Fig. 17. TT, EE, and BB angular power spectra of Madam (with 1 -minute baselines and without the noise filter) and MapCUMBA 
residual maps (black and light blue curves, respectively). Their RMS values were shown in Fig. [16] The maps represent the obser- 
vations of four LFI 30 GHz detectors (CMB, instrument noise, cooler fluctuations; asymmetric beams; and sample integration on). 
The red curve is for the difference of the Springtide (SPT) and MapCUMBA (MCB) residual maps. The blue curve is the same for 
Madam of 1 -minute baselines (MDM ( 1 min)) and MapCUMBA, and the green curve is for Madam of 1 .2 s baselines (MDM ( 1 .2 s)) 
and MapCUMBA. Madam used either 1.2 s or 1 -minute baselines, whereas Springtide used long 1 hour baselines. For comparison 
we also show the CMB spectrum of our simulations (gray curve), the spectrum of the Madam signal error map (including all sky 
emissions, asymmetric beams, and sample integration, dotted curve), and a theoretical CMB B-mode spectrum (dashed curve, see 
also Fig.[T2l. 

TABLE 3 



Statistics of signal error maps 11 



Code 




/ 






Q 






U 




MIN 


MAX 


RMS 


MIN 


MAX 


RMS 


MIN 


MAX 


RMS 


Madam (1 min) b 


-6.1 


3.5 


0.415 


-1.7 


1.5 


0.137 


-2.4 


1.6 


0.181 


Madam (1.2 s) b 


-49.4 


20.5 


0.890 


-7.8 


6.5 


0.307 


-11.3 


7.9 


0.421 


MapCUMBA 


-56.4 


21.5 


0.910 


-8.5 


7.0 


0.315 


-12.5 


8.9 


0.430 


ROMA 


-56.4 


21.5 


0.910 


-8.5 


7.0 


0.315 


-12.5 


8.9 


0.430 



11 We show here the statistics of the signal error maps (given in yuK) of different mapmaking codes. 
This table is for the case of asymmetric beams and sampling off. The corresponding maps are 
shown in Fig.[T8] They contain all sky emissions (CMB, dipole, and foregrounds). 

b We used Madam in two different configurations: with 1 -minute baselines and no noise filter 
(Madam (1 min)), and with 1.2 s baselines and with noise filter (Madam (1.2 s)). 



other Stokes parameters (Q and U) were larger as expected, but 
behaved similarly as a function of baseline length. 

The RMS does not provide a complete comparison of map 
quality. It is weighted toward the high-£ part of the spectrum, 



where in any realistic experiment we will be dominated by beam 
uncertainties and detector noise. Also, since all the errors in the 
map are folded into a single number, it tends to obscure the ori- 
gin of the errors. Fig. [17] shows typical angular power spectra 
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Fig. 18. Signal error maps of different mapmaking codes. The signal error is the difference between the noiseless output map and the 
binned noiseless map. These maps contain all sky emissions (CMB, dipole, and foregrounds) and they are in ecliptic coordinates. The 
beams were asymmetric and sample integration was off. The left column is for the temperature; the right column is for polarization 
magnitude. The statistics of these maps are given in Table 3. An example of the effects of beams and sample integration in the signal 
error maps is shown in Fig. [19] 



of the residual maps. In this plot the Madam (in this case with 
1 -minute baselines and without noise filter) and MapCUMBA 
residual spectra are nearly the same. The corresponding spec- 
tra of the other mapmaking codes would fall close to them too. 
To highlight the differences we show the spectra of three differ- 



ence maps between pairs of residual maps. Green curves show 
that the difference between Madam (with short baselines) and 
optimal residual maps is small. Comparison of the blue and red 
curves shows that destriping with 1 -hour baselines (Springtide) 
produces RMS residuals that are larger than those from the opti- 
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Fig. 19. Effects of beams and the sample integration in the signal error maps. We use Madam maps (1.2 s) here as examples. The 
maps contain all sky emissions and they are in ecliptic coordinates. The left-hand column is for the temperature and the right-hand 
column is for the polarization magnitude. The rows from the top correspond to symmetric beams & no sampling (SN), symmetric 
beams & sampling (SS), asymmetric beams & no sampling (AN), and asymmetric beams & sampling (AS). The statistics of these 
maps are given in Table 4. The third row maps are the same as the second row maps of Fig. [18] 



mal codes, and larger than those from destriping with 1 -minute 
baselines (Madam with 1 -minute baselines); the differences are 
confined to the high-/" part of the spectrum. At low i, the codes 
perform almost identically. 



To examine the signal errors of our mapmaking codes we 
made noiseless maps from the TOD of CMB, dipole, and 
foreground emissions. We subtracted the corresponding binned 
noiseless maps from noiseless output maps and obtained the 
signal error maps. We show some of our signal error maps in 
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TABLE 4 



Statistics of Madam signal error maps^ 



TOD 




/ 






Q 






U 




MIN 


MAX 


RMS 


MIN 


MAX 


RMS 


MIN 


MAX 


RMS 


Symm & no sampling 


-32.9 


14.3 


0.578 


-7.5 


6.1 


0.266 


-6.1 


5.0 


0.267 


Symm & sampling 


-30.2 


13.1 


0.609 


-7.3 


6.1 


0.272 


-5.8 


5.0 


0.272 


Asymm & no sampling 


-49.4 


20.5 


0.890 


-7.8 


6.5 


0.307 


-11.3 


7.9 


0.421 


Asymm & sampling 


-47.5 


21.9 


0.936 


-7.6 


6.5 


0.315 


-11.7 


8.1 


0.432 



a This table shows the effects that beams and sample integration have in the signal error maps. We show 
here the statistics of Madam (1.2 s) signal error maps (given in //K). The corresponding maps are shown 
in Fig. Q30 They contain all sky emissions (CMB, dipole, and foregrounds). The third line of this table 
is the same as the second line of Table 3. 




Fig. 15. Sample integration effects in the effective TT window 
functions. These window functions are the ratios of the TT angu- 
lar spectra of the binned noiseless CMB maps and the smoothed 
CMB input spectrum. The red curve is the same as the red curve 
of Fig. [14] The window functions of the symmetric beams (black 
and green curves) blow up at t ~ 600 due to pixelization error. 
The corresponding blow up of the other window functions oc- 
curs outside the scales of this plot. The thin black curve that 
tracks the green curve was computed from Eq. 0. The other 
black curve (that tracks the blue curve) is the window function 
prediction of our analytical model multiplied with Eq. 0. 



Figs. rT~8l and [T9l and their statistics in Tables 3 and 4. In Fig.fTSI 
we show signal error maps of a number of mapmaking meth- 
ods. These maps were made with asymmetric beams, with sam- 
pling off. Their statistics are given in Table 3. Optimal codes 
and Madam with short baselines (1.2 s) produce nearly the same 
signal error, which is stronger than the signal error of destripers 
with long baselines (represented here by another Madam map 
with 1 -minute baselines and no noise filter this time). For opti- 
mal and Madam (1.2 s) maps the signal error is more localized 
to the vicinity of the galaxy (which has strong signal gradients) 
than for long-baseline destriper maps. Fig. [17] shows that in the 
high-resolution 30 GHz maps the signal error is a small effect 
compared to the residual noise or the CMB signal that we used 
in this study. These results are well in line with the results of 



Fig. 16. RMS of the residual temperature maps for LFI 30 GHz 
observations. The maps were N^fe =512 full sky maps and they 
contained CMB, instrument noise, and the effects of the sorption 
cooler temperature fluctuations. Asymmetric beams were used, 
and sample integration was turned on. The y-axis units are ther- 
modynamic microkelvins. The angular power spectra of some of 
the residual maps are shown in Fig. [171 



our earlier studies (Poutanen et al. 2006, Ashdown et al. 2007b). 
In the bottom panel of Fig. [17] we compare the Madam signal 
error with the theoretical CMB B-mode spectrum (10 % tensor- 
to-scalar ratio). The plot shows that the magnitude of the sig- 
nal error is comparable to the magnitude of this B-mode signal. 
Signal error can therefore limit our possibilities to detect it. In a 
previous study we examined a number of techniques to decrease 
the signal error (Ashdown et al. 2007b). Because the detection 
of the B-mode signal was not a goal of this paper we did not 
investigate these methods for this data. 

Figure[19]and Table 4 show the effects of beams and sample 
integration on the signal error maps. We use the Madam maps 
here as examples. The figure and the table show that both switch- 
ing on the sample integration and switching from symmetric to 
asymmetric beams increase the signal error. This is because with 
asymmetric beams also the beam orientation affects the mea- 
sured signal. 
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Finally, we turn again to the issue of asymmetric beams, now 
in the presence of CMB and detector noise. Noise dominates the 
EE and BB spectra of our observed 30 GHz CMB maps. The TE 
spectrum, however, does not have a noise bias. Therefore, the 
effect of the T — > P cross-coupling can be detected in the noisy 
TE spectrum. We made a Madam map (with 1 -minute baselines 
and no noise filter) from the TOD of CMB+noise (asymmetric 
beams and no sample integration). In the TE spectrum of this 
map, the bias due to T — > P (arising from the beam mismatch) 
is clearly visible (see Fig. |20l ). We can expect that such a large 
bias will lead to errors in the cosmological studies (e.g., in cos- 
mological parameters), if not corrected. Fig. 1201 shows that our 
analytical correction method developed in AppendixlAlis able to 
restore the spectrum, at least on medium and large angular scales 
(at t < 450). 

5.3. Cooler fluctuations 

Temperature fluctuations of the Planck sorption cooler were de- 
scribed in Sect. 13.4.21 These fluctuations have an effect in the 
output signals (TOD) of the LFI 30 GHz radiometers. The typi- 
cal cooler TOD waveform was shown in Fig. [4] The RMS of this 
cooler signal is -35 fiK, which is about 1/38 of the RMS of the 
random uncorrelated instrument noise (white noise). 

It is of interest to bin the one-year cooler TODs of all four 
30 GHz detectors (identical TODs as described in Sect. 13.4.21 ) in 
a map, shown in Fig. [21] It looks similar to a map of correlated 
random noise (with faint stripes along the scan paths). Nutation 
of the satellite spin axis and the fluctuation of its spin rate (see 
Sect. 13. Il l randomize the regular cooler TOD signal when we 
project it in the sky. Therefore all map structures that these reg- 
ularities could produce are washed out. Because a pair of detec- 
tors sharing a horn see the same cooler signal, we might expect 
no cooler effect in the polarization maps. This is not, however, 
the case in reality. A small polarization signal arises because 
the polarization axes of the detectors are not exactly orthogonal 
within a horn (see Table 1). 

We made maps of CMB+noise+cooler and CMB+noise. We 
computed their difference to see the residuals of the cooler fluc- 
tuations that our mapmaking codes leave in their output maps. 
We computed the angular power spectra of these residual maps 
and show them in Fig. [22] Except for Springtide at small an- 
gular scales the residuals of the cooler fluctuations are smaller 
than signal error. We can therefore conclude, that, in these sim- 
ulations, the cooler effect is a tiny signal compared to the CMB 
itself, or to random instrument noise. 

The period of the cooler signal is -12 min (Fig. @). 
Springtide with its long (1 hr) baselines is not able to fit out 
any cooler power from the observed TOD. Therefore the cooler 
signal is not suppressed in the Springtide temperature map, 
but remains in the same level as in the binned cooler map. 
Springtide uses noise estimates that it computes from its ring 
maps. Although the TODs of this study have the same noise 
spectra the white noise levels that Springtide estimates for the 
rings will be different. The rings are combined using inverse 
noise variance weighting, and so when rings from two detec- 
tors in the same horn (which see the same cooler signal) are not 
weighted exactly the same, the effect in polarization is enhanced 
over the binned case which is generated using equal weights. We 
expect that Springtide would deliver EE and BB cooler spectra 
that are similar to the binned case if equal weighting would be 
used. 

Madam and the optimal codes, which operate with short 
baselines (optimal codes have one-sample-long effective base- 



Binned cooler T map 




Fig. 21. Temperature map binned from the cooler noise TOD. 
The map is in ecliptic coordinates and it contains observations 
of all four LFI 30 GHz detectors. Its statistics are (MAX, MIN, 
RMS) = (10.9, -11.7, 1.23) juK. We do not show the correspond- 
ing Q and U maps, but their RMSs were 0.002 and 0.007 //K, 
respectively. 

lines), are able to set their baselines to track well the cooler sig- 
nal and therefore the cooler signal is suppressed in their output 
maps. Our optimal codes and Madam assumed the same noise 
spectrum for all TODs. None of them tried to include the cooler 
spectrum in their noise filters. 

5.4. Pointing Errors 

We modified the Madam destriper in order to quantify how inac- 
curacy in the knowledge of satellite pointing affects the maps. To 
isolate this from other effects, we considered first the symmetric- 
beam, no-sampling case. However, to see the effect on polar- 
ization leakage, we needed a reference TOD which contained 
Stokes I only. This was available only for the asymmetric -beam 
no-sampling case (see Sec. 15.1.11 . 

We did not have available a good model for the expected 
pointing error, but to quantify the magnitude of the effect we 
considered two opposite schemes: uncorrelated pointing error 
with a random offset added to the pointing of each detector sam- 
ple and strongly correlated error for which the offset was kept 
constant for an hour at a time. The first model is of course unre- 
alistic, but serves to indicate to which direction the effects may 
change when the offset varies over shorter time scales. 

The pointing errors were drawn from a Gaussian distribution 
so that the RMS error for both the 6 and (f> directions was <x = 
0'.5. This translates into RMS offset of V2cr. 

The errors were generated independently for each horn; but 
the two detectors (for the two polarization directions) in the same 
horn share the same pointing and thus the same pointing error. 
Thus to a first approximation we would not expect the pointing 
error to contribute to the temperature-to-polarization leakage. 
However, since the two polarization directions were not aligned 
at exactly 90° to each other, there will be some effect. The differ- 
ent pointing errors of the two horns should contribute to E — > B 
leakage. 

The cr = 0'5 pointing error led to misplacing 11.5% of the 
samples into wrong pixels (Af s jd e = 512 pixels, about 7' across). 

In the limit of an infinitely dense, homogeneous hit distribu- 
tion, a simple model for the pointing errors can be devised. The 
errors increase the effective size of the pixels, leading to extra 
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Fig. 20. Effects of beam mismatch on the noisy TE spectra. Left hand panel: TE spectrum of a noisy CMB map (Madam map with 
1 minute baselines and no noise filter), black curve. To reduce its £ to t variation, it was filtered with a sliding average (A£ = 20). TE 
spectra of the binned noiseless CMB map (blue curve) and the smoothed input map (red curve) are shown too. They are the same 
curves as in the right-hand panel of Fig. [TT] The maps were made with asymmetric beams, and sample integration was off. Right 
hand panel: Noisy Madam TE spectrum (black curve), where the effects of beam mismatch have been corrected using the analytical 
method that we developed in Appendix. lAl 



smoothing of the output map, similar to smoothing with an ax- 
ially symmetric Gaussian beam with FWHM V8 ln2cr. In the 
power spectrum level this means an extra suppression factor of 
exp(-cr 2 ^+ 1)). 

We defined earlier the residual map as the difference between 
the output map and the binned noiseless map. Since pointing 
errors affect also the binned noiseless map, it is not enough to 
consider the residual map. Therefore we consider here the total 
error defined as the difference between the output map and the 
smoothed input map. 

Actual measurements typically include inhomogeneous dis- 
tributions of hits per pixel and within a pixel. Addition of ran- 
dom pointing errors can artificially smooth this distribution, re- 
distributing noise more evenly among the pixels. The expected 
white noise contribution to the RMS total error is proportional 
to V(l/Wobs)- Random pointing errors may decrease (1/A^ b s ). 
Thus random pointing errors may actually reduce the RMS total 
error in the map. However this does not happen with correlated 
pointing errors. 

Isolated from the redistribution of white noise are the effects 
on the actual signal. We divide the signal contribution to the total 
error into two parts: 1) the pixelisation error from uneven sam- 
pling of the pixel, defined as the difference (binned noiseless 
map - smoothed input map); and 2) signal error from destrip- 
ing, defined as the difference (destriped noiseless map - binned 
noiseless map). Pointing error affects both differences. The ef- 
fects of the pointing errors on these different error components 
are shown in Table 5 for temperature maps. Interestingly, un- 
corrected pointing errors also reduced the pixelization error; al- 
though the pixel was now sampled from a larger area, apparently 
the sampling was now more uniform. 

The total effect of the cr = 0.'5 pointing error was relatively 
small. In the uncorrected case, the RMS total error in the map 
decreased from 43.25 yi/K to 43.21 pK. In the correlated case it 
increased to 43.30//K. These effects are smaller than the dif- 
ference between the optimal codes and the destriper codes with 
long baselines (i.e., Madam without noise filter, Springtide). 



We then studied the angular power spectra of noiseless maps 
binned with additional pointing errors. Should the naive model 
derived from an infinitely dense homogeneous hit distribution 
hold, the only difference between spectra from different cases 
would be the added smoothing. We found this to be the case 
approximately when binning only the CMB signal. Adding dif- 
fuse foregrounds and point sources broke the approximation to 
a large part, presumably due to sharper structure and lack of sta- 
tistical isotropy in these foreground components. The decreasing 
trend from extra smoothing remained still visible. In all cases 
the random pointing error model corresponded to the smoothing 
approximation considerably better. The spectra are displayed in 
Fig. [23] Pointing errors contribute also to the aliasing effect {(- 
mode coupling due to uneven distribution of hits), which adds 
power to high £, but the effect is negligible in comparison to the 
residual noise. 

For analyzing the effect on polarization maps and power 
spectra we repeated the study using the asymmetric beam TOD, 
since for that case a reference TOD containing Stokes I only was 
available. See Table 6. For correlated pointing errors the increase 
of the residual noise error is larger by almost an order of magni- 
tude for polarization than for temperature (Table 5). This effect 
is in line with the estimate based on white noise. The effect is 
probably due to a relatively small number of pixels where the 
sampling of polarization directions is not very good, and thus 
the noise effects get magnified for polarization. The other effects 
on the map total RMS error are of the same magnitude for polar- 
ization maps as for temperature maps. 

The T — > P and E — > B leakage was studied the same way 
as in Sect. 15.1.11 The effect of pointing error on this leakage is 
shown in Fig. [24] The effect on T — > E leakage is negligible, and 
the largest effect is on the E — > B leakage, as anticipated. This 
should be compared to the magnitude of the kind of B mode 
we might expect to detect with Planck. In the right-hand panel 
of Fig. [24] the leakage is compared to the theoretical fi-mode 
spectrum that we introduced earlier (see Fig.[T2l. It can be seen 
that the leakage to B is small compared to this signal for £ < 
300. On the other hand, the pointing error effect on this leakage 
becomes small for £ > 200. 
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Fig. 22. Effect of sorption cooler temperature fluctuations on the CMB maps. These plots show the TT, EE, and BB angu- 
lar power spectra of residual maps of the cooler effect. The residual maps were computed as differences between the maps of 
CMB+noise+cooler and CMB+noise. They represent one year of observations of four LFI 30 GHz detectors. The beams were 
asymmetric and sample integration was on. The light blue curves ("Binned") show the spectrum of the binned cooler map, whose 
temperature component is in Fig. [21] Note that the MapCUMBA, ROMA, and Madam curves are nearly on top of each other. For 
comparison we also show the spectrum of the Springtide CMB+noise+cooler map and the spectrum of the Madam signal error map 
(including all sky emissions, asymmetric beams, and sample integration, dotted curve, same curve as in Fig.fTTIi. 



TABLE 5 

Effect of satellite pointing errors on total temperature error maps 2 









Error components 


(T 


Model 


Total error White noise estimate 


Residual noise b Signal error c Pixelization error d 


0' 

0'.5 
0'.5 


uncorr. 
corr. 


43.2560 42.5449 
-0.0432 -0.0279 
+0.0482 +0.0222 


42.6536 0.5601 7.1663 
-0.0318 +0.0023 -0.0611 
+0.0166 +0.0383 +0.2047 



a We show here how the total RMS error (given in //K) in the map changed by adding pointing errors. This 
table is for the case of symmetric beams. The first row shows the case without pointing errors. The two other 
rows show how the map errors changed when adding pointing errors according to the two models considered 
in Sec. 15.41 The total map error is divided into three parts: residual noise, signal error, and pixelization error, 
and the effect of pointing error to each part is shown. 

b Residual noise is the power in the difference (destriped map - destriped noiseless map). It is dominated by 
binned white noise. 

c Signal error is the destriping error caused by the signal (destriped noiseless map - binned noiseless map). 
d Pixelization error is caused by insufficient sampling of the pixel temperature. We measure it from (binned 
noiseless map - smoothed input map). 



Based on these results, pointing errors of < 0'5 do not appear 
as a major concern for the 30 GHz channel. 
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TABLE 6 

Effect of satellite pointing errors on total polarization error maps" 











Error components 




Model 


Total error 


White noise estimate 


Residual noise b Signal error Pixelization error 


0' 

0'.5 
0'.5 


uncorr. 
corr. 


85.4361 
-0.0444 
+0.1591 


85.0429 
-0.0543 
+0.1428 


85.3790 0.3644 3.4872 
-0.0398 +0.0024 -0.0949 
+0.1364 +0.0439 +0.5112 



a We show here how the total RMS error in the polarization maps changed by adding pointing errors. 
Polarization error amplitude is the square root of the square sum for the two polarization amplitudes: 

\err po i\ = ^jerr 2 Q + err 2 ,. 

b See Table 5 for the descriptions of the columns. 
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Fig. 23. TT power spectra of the 0'5 pointing error cases divided by the corresponding spectra without pointing errors. All spectra 
are for the binned noiseless maps (fg = foreground, ps = point source). These plots show that for the CMB only case the pointing 
errors are well approximated by an effective pixel window (due to larger effective pixel size). However the simple model fails in the 
presence of secondary signals. 



6. Computing resources 

Timing tests of codes were run on the NERSCpl computer 
Jacquard, an Opteron cluster running Linux, chosen because it 
is similar in architecture to the present and future machines in 
the Planck Data Processing Centers (DPCs). Results are shown 
in Table 7 and in Fig. [25] Different codes have varying demands 
for processing power, memory, network interconnect, and disk 
I/O. Different computers will provide these resources with vary- 
ing performance, which can make a significant difference in the 
runtime of a particular code. For this reason, the Jacquard tim- 
ing tests should be viewed as a useful metric, but not an absolute 
prediction of performance. Note however that many of the codes 
perform similar tasks (for example, each must read the TOD), 
so a change in the performance between machines (in this exam- 



12 National Energy Research Scientific Computing Center, 
http://www.nersc.gov 



pie, in disk I/O) will have less impact on the relative performance 
measured between codes. In the real-world application of these 
codes, they will be tuned to the particular hardware available. 

In early 2009, the LFI DPC will have at least 128 processors 
with 1 GB per processor (128 GB total memory), with a goal of 
256 processors and 2 GB per processor (5 12 GB total memory), 
plus 10 TB of disk space. The HFI DPC expects to have 250 pro- 
cessors and 8 GB per processor (2000 GB total memory), plus 
200 TB of disk space. 

To compare the fast mapmaking capabilities of our codes, we 
used the data of Table 7 to calculate the rate that 30 GHz maps 
could be produced in 24 hours, assuming that 50% of the pro- 
cessors and the memory of the LFI DPC goal computer would 
be available. Results are shown in Table 8. The optimal codes 
(MADmap, MapCUMBA, ROMA) could produce a couple of 
hundred maps per day, compared to between several hundred 
and a couple of thousand for the destripers (Springtide, Madam). 
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Fig. 24. Polarization angular power spectra of the maps made from the noiseless temperature-only TOD and the difference maps 
between these maps and the ones made from the noiseless TOD with polarization, showing the E-mode polarization signal, and the 
T — > E, T — > B and the E — > B leakage effects. For comparison, we show a theoretical B mode spectrum corresponding to a 10% 
tensor-to-scalar ratio and including lensing from E. 

TABLE 7 

Resource Requirements 3 

Memory [GB] Time 



CODE #Proc. Per Processor Total Wall Clock [s] Total [CPU-hr] I/O [CPU-hr] 



Springtide (A b ) 16 1.0 16.1 293 1.3 

Springtide (B) 16 0.8 13 336 1.5 

Madam (1-hr baseline, A) 16 0.6 9.4 176 0.8 0.27 

Madam (1-hr baseline, B) 16 0.2 2.7 318 1.4 0.64 

Madam (1.2s baseline, A) 16 2.2 34.8 362 1.6 0.23 

Madam (1.2s baseline, B) 16 0.5 8.2 1167 5.9 0.29 

MADmap (A) 96 1.7 166 427 11.4 0.25 

MADmap(B) 32 1.1 35.2 1953 17.4 0.71 

MapCUMBA 64 1.9 122 1441 25.6 9.73 

ROMA 64 2.1 134.4 1590 28.3 15.9 



a Results are for runs on the NERSC Jacquard computer. Typical run-to-run variations in total CPU-hours (central processing 

unit) are 10%, determined from Madam runs. In all tests, N s a c = 512 (/, Q, U) maps were made from 12 months of simulated 

observations from four LFI 30 GHz detectors. 
b Several codes can trade processing cost for memory usage. For a single code, we mark the timing run optimized for CPU by 

(A) and the run optimized for memory by (B). 
c In Springtide, reading the TOD from disk occurs simultaneously with the compression of the TOD to ring-maps. This makes a 

separate evaluation of the I/O difficult. 



For the LFI and HFI DPCs, the most challenging mapmaking 
problem is to make a map from 14 months of observations (the 
nominal Planck mission time) using 12 detectors (the maximum 
for any channel) at 70 and 217 GHz. We have estimated the total 
memory and disk space required by the mapmaking codes in 
these extreme cases. 

For the LFI 70 GHz channel, a naive scaling up of the size of 
the TOD from the 30 GHz timing tests (considering number of 
detectors, sampling rate, and mission duration) yields total mem- 
ory requirements 8.3 times those shown in Table 7. The Madam 
team made more detailed estimates based on counting the size of 
allocated arrays. For the code running in standard configuration 



with 1.2-s baselines, the estimate yields a total memory require- 
ment of ~ 230 GB, compared to the naive scaling estimate of 
270 GB. In split-mode, where the data are first destriped with a 
1.2-s baseline in three groups of four detectors, then combined 
in a second destriping using 1-hr baselines, the memory require- 
ment is 100 GB. For Madam with 1-hr baselines, the memory 
estimate is 30 GB. In split-mode and for 1-hr baselines, the esti- 
mate is less reliable because of uncertainties in the performance 
of Madam's data compression system. For a code using the full 
detector pointing, the disk space required to store the TOD is 
0.87 TB. For the compressed satellite pointing sampled at 1 Hz, 
the requirement drops to 0. 12 TB. 
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TABLE 8 



LFI 30 GHZ MAPMAKING THROUGHPUT" 



Code 


Max. Simultaneous Maps 


Maps/day 


Springtide (A) 


8 


2359 


Springtide (B) 


8 


2057 


Madam (1-hr baseline, A) . . . . 


8 


3927 


Madam (1-hr baseline, B) . . . . 


8 


2173 


Madam (1.2 s baseline, A) . . . . 


4 


954 


Madam (1.2 s baseline, B) . . . . 


8 


592 


MapCUMBA 


2 


119 


MADmap (A) 


1 


269 


MADmap(B) 


4 


176 


MapCUMBA 


2 


119 


ROMA 


1 


108 



a Assumptions: 1) 128 processors with total 256 GB RAM are allocated to map- 
making (~50% of the LFI DPC goal computer); 2) code performance is taken 
from Table 7; 3) computer resources are allocated in 16-processor nodes. 

b Obtained by dividing the number of seconds in a day by the cpu usage in Table 7, 
taking into account the number of utilized nodes, processors per node, and num- 
ber of simultaneous maps. 
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Fig. 25. Codes with longer baselines can run with fewer re- 
sources. Individual codes can be tuned to be compact in memory 
or to run quickly. 



For the HFI 217 GHz channel, the naive scaling yields mem- 
ory requirements 21.5 times those shown in Table 7. The Madam 
detailed estimate yields 1070 GB for 1.2-s baselines, 460 GB for 
1.2-s baselines in split mode, and 100GB for 1-hr baselines. 
Using full pointing, the TOD disk requirement is 2.25 TB, but 
drops to 0.32 TB for compressed pointing. 



7. Conclusions 

In this paper we have presented results from a simulation that 
we designed to determine how mapmaking codes handle four as- 
pects of real Planck data: non-axially-symmetric beams, detec- 
tor sample integration, temperature fluctuations induced by the 
sorption cooler, and pointing errors. Five mapmaking codes of 
two basic types (destriping and optimal) were used in this study. 



We generated one year long streams of observations representing 
four LFI 30 GHz detectors. This simulation round (Trieste sim- 
ulation) was number four in a series of mapmaking comparison 
studies that the Planck CTP group has undertaken. 

None of our mapmaking codes made an attempt to de- 
convolve the beam effects from its output map. Therefore the 
smoothing effects of the beam and sample integration showed 
up similarly in all our maps. This is, however, a complicated 
smoothing, because every map pixel has its own effective beam. 
We also made a thorough examination of the temperature and 
polarization cross-couplings that arose from the mismatch of the 
beams. These effects were also similar in all our maps. 

Our conclusions are: 



Our studies showed that the temperature to polarization 
cross-coupling of CMB signal caused a detectable bias in 
the TE spectrum of the CMB map. In the EE spectrum the 
effect was small compared to the residual noise of the map. 
The E — > B cross-coupling of CMB produces a spurious sig- 
nal whose magnitude in intermediate and large multipoles 
(at i > 300 in this case) seems to exceed the magnitude of 
the CMB fi-mode signal that we might expect to detect with 
LFI 30 GHz detectors. However, in this range of multipoles 
the detection of the B-mode is difficult due to map noise any- 
way, and at low multipoles where the signal-to-noise ratio is 
higher the effect is small. Our study shows that the spectrum 
biases that T — > P and E — > B cross-couplings cause may 
lead to errors in the cosmological studies, if not corrected. 
We show a method in this paper that can correct part of these 
effects in the spectrum domain. Map-domain methods that 
address the beam convolution issues properly have been pro- 
posed by other authors (see Sect. 15. 1.3] ). 
Signal error (error that mapmaking couples to the output map 
from the small-scale signal structure) is the part of the out- 
put map where beam effects appear differently in the maps 
of different mapmaking codes. In optimal and short baseline 
Madam maps this error is larger than in the destriper maps of 
long baselines (Madam and Springtide). Signal error is, how- 
ever, a small effect in a high-resolution map if we compare it 
to the residual noise of the map. Techniques to decrease the 
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signal error were thoroughly examined in our earlier study 
(Ashdown et al. 12007b) . 
3. Based on the results of this study, cooler noise or pointing 
errors of < 0'5 do not appear as a major concern for the 
30 GHz channel. 

Five mapmaking codes (two destripers and three optimal 
codes) have been developed and compared in four Planck CTP 
simulation rounds (Cambridge, Helsinki, Paris, and Trieste sim- 
ulations). We give the following summary of these studies: 

1. At the end of a long process, the only essential difference 
found so far between codes that affects accuracy — as as- 
sessed by the RMS residual between an input binned map 
and the output map — is baseline length, shown in Fig. [16] 
This assumes that all maps have been made at sufficiently 
high resolution that the effects of sub-pixel structure are 
small (see Ashdown et al. 2007b). (Optimal codes are some- 
what more affected [in an RMS residual sense] by sub-pixel 
structure than are destripers, but this is not an issue at suffi- 
ciently high map resolution.) However, the fundamental dif- 
ference between destriping and optimal codes is in their as- 
sumptions about, and handling of, noise. Some of the future 
tests listed in Sect. 17. ll will probe this difference more inten- 
sively than ones performed to date. 

2. All optimal codes give essentially indistinguishable results. 
Madam, a destriper with adjustable baseline length, gives the 
same result as the optimal codes when the baseline is set 
short enough. Reliable noise estimation is required for these 
codes. 

3. Resource requirements for the codes vary by more than an 
order of magnitude in the important categories (memory per 
processor, total memory, elapsed time, and total CPU-hours), 
as shown in Table 7 and Fig. [25] Codes set up to deal with 
pointing information in multiple forms have a practical ad- 
vantage, in that the resources required in memory, I/O speed, 
and processor speed can be matched to the resources avail- 
able (memory, I/O speed, and CPU time). 

4. It seems inevitable that Planck will need both a fast and "re- 
source light" mapmaking code for everyday use, and an "ul- 
timate accuracy" but "resource heavy" code, for final data 
products. 

7. 1 . Future Tests 

Our final round of simulations was quite realistic compared to 
earlier rounds, nevertheless there are still many instrumental and 
mission effects, both subtle and blatant, that have not been in- 
cluded, some of which are listed below. Further tests of map- 
making must be performed. 

- Gaps in TOD. 

- High-pass filtering of TOD (e.g., to reduce low frequency 
noise before mapmaking, to deconvolve the bolometer fre- 
quency response in HFI). High-pass filtering will make the 
noise matrix C„ of Eq. ([3]) non-diagonal. 

- Leakage of galactic temperature signal to the polarization 
maps due to the bandpass mismatch of the LFI radiometers. 

- Beam sidelobes. 

- Cross-polarization leakage of the detectors, and polarization 
angle errors. 

- Other forms of random noise (e.g., "popcorn noise"). 

- Features in the noise power spectrum (e.g., microphonics, 
frequency spikes). Note that some of these features may 
make the noise matrix C„ of Eq. ([3]) non-diagonal. 



- Changes in the baseline signal level (e.g., from telescope 
emission, background loading of bolometers). 

- Realistic pointing errors (errors in the pointing reconstruc- 
tion) 

The scope of this paper was restricted to the mapmaking 
from the LFI 30 GHz data. The Planck mapmaking codes need, 
however, to be able to make maps from both LFI and HFI ob- 
servations. Therefore the list includes mapmaking tests that are 
relevant for LFI, HFI or both. 

Future tests should be run on both a destriping code and an 
optimal code. As a practical matter, it will be efficient to run 
the less resource intensive destriping code first to wring out any 
problems in setting up the test. It is no longer necessary to run 
tests on multiple implementations of destriping or optimal codes, 
which have been shown to give the same answer. 

Generally speaking, destriping algorithms assume that the 
noise power spectrum is white at high frequency. The shape of 
the low frequency power spectrum is not important, so long as 
the low frequency noise may be fit by a series of offsets (or 
baselines). The duration of these offsets varies from 1 second to 
1 hour in the codes we have tested here, with a trade-off between 
accuracy and resource usage. For baseline durations longer than 
1 minute, destripers need no information about the noise power 
spectrum. By contrast, the optimal codes are not restricted to 
white noise at high frequency, but must be informed of the shape 
of the noise power spectrum. 

Deconvolution of a bolometer time constant exercises a dif- 
ference in these approaches. Bolometers have a response lag to a 
change in the sky similar to the time for the beam to cross a fixed 
point. This effect smears the sky image in the direction of the 
scan. A time stream deconvolution filter, applied prior to map- 
making, can repair this effect on the sky signal, but at the cost 
of correlating the noise at small scales. We have tested the effect 
of this with simulated HFI 100 GHz data, where the bolometer 
time constants are the longest, and where the effect is the most 
severe (Gorski et al. 12008) . 
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Appendix A: Effects of beams in the angular power 
spectrum of a CMB map, an analytic model 

In this Appendix we try to explain the mechanisms by which 
beam mismatches generate cross-couplings between the Stokes 
/, Q, and U of a CMB map. We designed a simple analytic model 
to answer this question. We start by making some simplifying 
assumptions: 



- Sky pixels are observed same number of times by all four 
LFI 30 GHz detectors. 

- All hits fall in the centers of the pixels. 

- Detector polarization angles are uniformly spaced in 180° 
(spacing is 45° for the four 30 GHz detectors). 

- Detector noise is white, Gaussian, and every detector has the 
same noise RMS. 

- We consider only the main beams. Far sidelobes are ignored. 

- The widths of the main beams are so small that we can ap- 
ply the small-scale limit and compute the convolution of the 
beam and sky in the plane wave Fourier space (flat sky ap- 
proximation). 

- Our asymmetric main beams have elliptical Gaussian re- 
sponses. 

The assumptions lead to the following consequences: 

- Ellipses are invariant to 180° rotations, therefore the scan- 
ning direction does not matter. The results of our model are 
the same in north-to-south and south-to-north scans. 

- In all parts of the sky a detector beam has the same orienta- 
tion with respect to the local meridian. 

- The model does not account for the effects of the pixel sam- 
pling (pixel window function and pixelization error). 

We expect that in spite of the simplifying assumptions, our 
model gives a good description of the observations in most parts 
of the sky. The accuracy of the model should be good in those 
parts of the sky that are scanned approximately along meridi- 
ans (ecliptic equator areas). The accuracy of our model may be 
worse in the ecliptic pole areas, where there is more spread of 
the scanning directions. In reality, the polar regions with multi- 
ple scanning directions are a relatively small fraction of the total 
sky (see Sect. |5]and the third row of Fig. |5). 

The Stokes parameters Q and U at a point in the sky are 
defined in a reference coordinate system (e^e^, n), where the 
unit vector e<j is along the north-south meridian (increasing 6), 
e<£ is along the increasing (p, and n points to the sky (Gorski et 
al. |2005bl >. In this reference coordinate system the polarization 
angle of a detector is the angle from the local north-south merid- 
ian to its polarization sensitive direction (anti-clockwise rotation 
in the reference coordinate system). The polarization angles of 
the four detectors (LFI-27a, LFI-27b, LFI-28a, LFI-28bl3 are 
ij/, ij/+90°, i//+45°, and <//+135°. In the north-south scanning the 
value of i/> is the same as the angle from the scanning direction 
to the U-axis of LFI-27a (see Fig. |2j. This angle is 67?5 (see the 
caption of the same figure). For the opposite scanning direction 
(A is increased by 180°. 

The observations from a sky pixel p can be given as 



y(n p ) = 



(yi(n P )) 


1 






y3(n p ) 




, y4(n p ) , 


\ 




/ 




V 



h+Qx cos(2t/0 + Ui sin(2i/r) N 
h - Qi cos(2t/>) - U2 sin(2i/0 
h - Q-i cos(2i/0 + U 3 sin(2i/0 
I4 + Q4 cos(2i/i) - U4 sin(2i/f) 

h +e- a * ■Zi/2 + e +i2 'f' -Z*/2 
h - e-W ■ Z2/2 - e +il * ■ Z*/2 
h - ie~ i2 * • Z 3 /2 + ie+°* ■ Z*/2 
I A + ie~ i2 * ■ Z4/2 - ie +i2 fZ* 12 ) 



(A.l) 



Here indices (1,2,3,4) refer to the detectors (LFI-27a, LFI-27b, 
LFI-28a, LFI-28b), X* is the complex conjugate of X, Z = Q+iU 
is the complex polarization field, and Z* = Q — iU. In our model 



- Scanning is from one ecliptic pole to another along the 
meridians. 



13 Detectors LFI-27a and LFI-27b share a horn and the detectors LFI- 
28a and LFI-28b share the other horn. 
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the detectors measure different values of the Stokes parameters 
from the same point of the sky, because their beams are differ- 
ent. The observations stay the same in both scanning directions, 
because the beams have 180° rotational symmetry and a change 
of iff by 180° does not change anything in Eq. (IA.U . 
Let us define a "pointing matrix" A 



A = 



1 



e 



1 -ie-^ ie +i2 * 
1 ie-^ -ie +i2 * 



(A.2) 



and a Stokes triplet s ; , of the pixel p 



s„ = 



Kn P ) 
Z(n p ) 
Z*(n p ) , 



(A.3) 



The unit vector n p points to the center of the pixel. 

For uniform white Gaussian noise, the Stokes map ( s p ) can 
be recovered from the detector observations as (Tegmark l 19971 1 



s p = (aIa) -1 A t y(n p ). 



(A.4) 



Here X' is the hermitian conjugate of matrix X. Matrix A' A is 
the N bs matrix of the pixel p (see the footnote of Sect. [5). The 
assumptions that we made in the beginning of this Appendix lead 
to diagonal Nobs (A 1 " A = diag[4, 1, 1]). 

We can now solve I(n p ) and Z(n p ) from Eq. (IA.4t . 



/(n P ) 



(7i + h + h + h) c 



-Hi// 



[(Z, -Z 2 )-i(Z 3 -Z 4 )] 



e +a f [(Z*-Z*) + i(Z;-Z*)] 



(A.5) 



Z(n„) 



(Zi + Z 2 + Z 3 + Z 4 ) | e + ' 2 ^ [(/i - h) + i(h - h)} 



[(Z* +z*)-(z; +Z*)] 



(A.6) 



The first terms on the right hand sides (of Eqs. ( IA.5b and 
dA.6l l) provide the temperature and polarization signals that we 
want. The other terms of Eq. (IA.5t represent the cross-coupling 
of polarization signal to the temperature. This cross-coupling 
can usually be ignored because the coupling coefficient is typ- 
ically small and the polarization signal is weaker than the tem- 
perature signal. 

The second and the third terms of Eq. jA.6i correspond to 
T — > P and spin-flip coupling errors (Hu et al. 2003 ). The sec- 
ond term shows that the temperature signal can pollute the po- 
larization signal if the beams of the two detectors sharing a horn 
do not match. Beams can be different from horn to horn, but 
no leakage of temperature to polarization occurs as long as the 
beams of a horn are identical. 

The spin-flip term (the third term of Eq. ( 1A.61 )) produces 
cross-coupling between E- and B-mode polarizations. We show 
it later in this Appendix. In our model spin-flip coupling occurs if 
there is a horn-to-horn mismatch between the total beams. Here 
the total beam means the sum of the responses of the two beams 
of a horn. 



Convolution of the detector beam and the sky pro- 
duces the observed temperature and polarization fields 
(Ii,Z\, I 2 ,Z 2 , 1 3 ,Zj,, I 4 ,Z 4 ). We assumed that this convolution 
can be computed as a multiplication between the correspond- 
ing flat-sky Fourier-domain quantities. This is a good approx- 
imation in a small flat patch around the point of observation. 
We can therefore write 7i(k) = Bi(k)/(k), Zj(k) = B[(k)Z(k), 
7 2 (k) = fi 2 (k)/(k), Z 2 (k) = B 2 (k)Z(k) and so on. Here k 
is the 2-dimensional plane wave vector and (Z?i(k), B 2 (k), 
Bj,(k),B 4 (k)) are the Fourier-domain representations of the (flat- 
sky) responses of the (LFI-27a, LFI-27b, LFI-28a, LFI-28b) 
main beams. Because we assumed ideal elliptic Gaussian main 
beams, their functional forms are well known and they are real- 
valued for all widths, ellipticities and orientations (Fosalba et 
al. [2002). The Stokes parameters of the sky are given by /(k) 
andZ(k). 

We can now write the relation between the Stokes parameters 
of the map and those of the sky (for Fourier-domain quantities): 



.(A.7) 



Here B s , Z?d, and B c are Fourier-domain quantities and they are 
defined in terms of the beam responses 
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B s (k) = Bj(k) + B 2 (k) + B 3 (k) + B 4 (k) 
B d (k) = B,(k) + B 2 (k) - B 3 (k) - B 4 (k) 
B c (k) = (Bi(k) - B 2 (k)) + i(B 3 (k) - B 4 (k)) . 
B* c (k) = (Bi(k) - B 2 (k)) - i (B 3 (k) - B 4 (k)) . 



(A.8) 
(A.9) 
(A. 10) 
(A. 11) 



The 3x3 matrix of the right-hand side of Eq. (IA.7t gives the 
mapping from the Stokes parameters of the sky to the Stokes 
parameters of our map. We can therefore call it the Mueller ma- 
trix. We use the symbol M(k) for it here. In a real experiment 
each pixel of the map has its own Mueller matrix, because de- 
tector hit counts, sampling of pixel area, and beam orientations 
are different for different pixels. Our model, however, leads to 
one Mueller matrix that applies in all pixels of the sky. 

Next we make the connection between the polarization field 
(Z(k)) and the fields of the E- and B-mode polarization (E{k) 
and B(k)) (Zaldarriaga & Seljak 1997). We include the tempera- 
ture anisotropy field (/(k)) in these equations. 
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(A. 12) 



Here the angle 0t is defined through k = (k x ,k y J = 
k (cos(0k), sin(0t)) and k = \k\ is the magnitude of the wave vec- 
tor. The definition of the 3x3 matrix R(0k) is evident from the 
equation. 

We assume that the sky emission (like the CMB) is statisti- 
cally isotropic. Therefore the ensemble mean of its power spec- 
trum does not depend on the direction of k and it can be com- 
puted in a 3 x 3 matrix form 



C(k) = 



(C TT (k) C TE (k) C TB (k) 
C TE (k) C EE (k) C EB (k) 
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.(A.13) 
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Here (X) is the ensemble mean of X. We can now compute the 
ensemble mean of the power spectrum matrix of the observed 
map as 



/(k) 
E(k) 
B(k) 



( too I 

E(k) 
B(k)J 



(A. 14) 



Eqs. ( IA.7l i and (IA. 12b give a relation between (/, E, B) of the 
sky and (/, E, B) of our map. With the help of this relation, we 
can write an equation that gives us the power spectrum of our 
map if we know the power spectrum of the sky. 



(C(k)) = ^J R'Vk) • M(k) ■ R(0 k ) ■ C(k) ■ 

■R(0 k ) ■ M + (k) ■ R 1 "^) • d<p k . (A. 15) 

The 3x3 power spectrum matrix is symmetric. It therefore 
contains 6 different component spectra C XY (k). We can arrange 
the power spectra of the input sky in a 6-element column vector 

( C TT (k) \ 
C TE (k) 
C TB (k) 
C EE (k) 
C EB (k) 
C BB (k) 



c(k) = 



(A. 16) 



Similarly we can build a 6-element column vector c(k) of the 
power spectra of the observed map. Because the power spectrum 
of the input sky does not depend on the direction of the wave 
vector, Eq. ( IA. 15b allows us to construct a 6 x 6 window matrix 
that, when applied to the spectrum c(k) of the input sky, gives 
the spectrum c(k) of the observed map. 

We can use the above results to predict the angular power 
spectrum (spherical harmonic domain) of our map, if we know 
the angular power spectrum of the sky. The predicted spectrum 
will include the effects of the beams, but it will not include the 
effects from the sampling of the pixel areas (pixel window func- 
tion and pixelization error). The steps to compute the angular 
power spectrum prediction are the following. 

1. Assume that we know the angular power spectrum multi- 
poles of the input sky. Arrange them in a 6-element column 
vector at every multipole t (as in Eq. ( lA.16b ). 

2. The magnitude k of the Fourier wave vector is the contin- 
uous limit of the integer t labeling the spherical harmonics 
multipoles a( m . 

3. For a given multipole (, set k — I in the Mueller matrix M(k). 

4. For the same multipole, step the angle 0k from to 2n and 
compute the 6x6 window matrix from the integral of Eq. 

(TaTT3T> . 

5. Apply the window matrix to the angular power spectrum 
vector of the input sky. The result is the angular power spec- 
trum vector of the observed map (of the same multipole €). 

6. Repeat the previous steps to all multipoles of interest. 

Temperature to polarization cross-coupling is an important 
effect of the beam mismatch. We can use the above procedure to 
predict this effect too. Instead of using the full 6-element input 
spectrum vector, we use a vector all of whose elements except 
C TT are zero, and carry out the prediction steps as before. 

We can invert the 6x6 window matrix and compute a cor- 
rected angular power spectrum from the spectrum of the ob- 
served map. In the corrected spectrum the beam effects have 



been deconvolved out and it is an approximation of the spectrum 
of the input sky. 

We made a software code that implements the above predic- 
tion and correction steps. Instead of building 6-element power 
spectrum vectors, we assigned the full 3x3 power spectrum 
matrix column-wise in the vector c(k). This leads to 9-element 
power spectrum vectors and 9x9 window matrices. They have 
unnecessary redundancy, but this approach makes the implemen- 
tation simple and straightforward. The extra computational cost 
of larger vectors and matrices is unimportant. 

Using Eq. ( IA. 15b we can compute the 9x9 window ma- 
trix (see Appendix A of Hamimeche & Lewis 2008 for relevant 
equations of matrix vectorization) 
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Here S = S (k) = R'l 
Eq. dA~T5l >) and S u 



6k) • M(k) ■ R(0 k ) is a 3 x 3 matrix (from 
= Sij(k) are its elements. The predicted 
spectrum of the observed map (from the spectrum of the input 
sky) becomes now 



c(k) = W(Jt) ■ c(k). 



(A. 18) 



We further convolved the predicted spectra with HEALPix pixel 
window function to model the smoothing due to sampling of the 
pixel area. 

The prediction equation (Eq. (IA.18I )) is simple to invert for 
the correction. Before applying the inverse W{k) to the spec- 
trum of the observed map, we deconvolved the observed spec- 
trum with the HEALPix pixel window function. 

We used this code in Sect. |5]of this paper and computed pre- 
dictions from the input spectrum and corrections from the map 
spectrum. These predictions we used to explain the beam mis- 
match effects and beam window functions that we detected in 
the spectra of our CMB maps. The correction capability of our 
model is demonstrated in Fig. [20] 

A. 1 . Spin-flip coupling 

If the (2,3) and (3,2) elements of the Mueller matrix are non-zero 
(see Eq. ( IA.7b ), Z* will mix with Z and vice versa. This is the 
spin-flip coupling (Hu et al. 2003). It generates a Z?-mode in the 
map even if the input sky has no B-mode power in it (like the 
CMB sky of this study). We can describe this situation with the 
following general expression (for the Fourier-domain quantities) 



Z 
Z* 



a b 
b* a 



Z 
Z* 



(A. 19) 



We assume that a is real as the (2,2) and (3,3) elements of the 
Mueller matrix (see Eq. jA.lj )- The quantity b and its complex 
conjugate represent the (2,3) and (3,2) elements of the Mueller 
matrix. If we assume that the input sky contains no B mode 
power (ie. Z and Z* arise from E mode only), we can write for 
the E and B modes of the mixed fields Z and Z* (use Eq. (1A.12I ) 
for the relation between the E and B modes and the complex 
polarization fields) 
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This equation can be put in the following form 

fc*g+<40„ + be -i4<l, k \ 

i(b*e +iA ^-be- iA ^)) E - (A ' 21) 

The equation shows that spin-flip coupling (non-zero b) creates 
B mode in the mixed field. In addition, the spin-flip coupling 
influences the E mode itself in the sense that E of b = and E 
of b + are different. 




